MATLAB Examples

# Horizontal line

This example presents the generation of a turbulent wind field along a horizontal line, assimilated to a bridge deck. The flow is assumed homogeneous and stationary along the deck.

## Run

The function windSim.m is called, with a single input. This input is a
text file named "INPUT.txt". In this example, this file is located in
the same folder as windSim.
clearvars;close all;clc; rng(1); % to ensure reproducibility of the example. filename = 'INPUT.txt'; [u,v,w,t,nodes] = windSim(filename); 
Elapsed time is 17.516536 seconds. 

## Time series overview

Visualization of the spatial correlation of the time series
figure subplot(311) box on axis tight plot(t(1:600),u(1:3,1:600)+nodes.U(1:3)*ones(1,600)) xlabel('time (s) '); ylabel('u (m/s)') set(gcf,'color','w') legend('sample 1','sample 2','sample 3'); subplot(312) box on axis tight plot(t(1:600),v(1:3,1:600)) xlabel('time (s) '); ylabel('v (m/s)') set(gcf,'color','w') subplot(313) box on axis tight plot(t(1:600),w(1:3,1:600)) xlabel('time (s) '); ylabel('w (m/s)') set(gcf,'color','w') 

## Turbulence length scales

Calculated from the auto-correlation functions.
dt = median(diff(t)); fs = 1/dt; tmax = t(end); f0 = 1/tmax; [Nsamples,N] = size(u); Lu = zeros(1,Nsamples); Lv = zeros(1,Nsamples); Lw = zeros(1,Nsamples); for ii=1:Nsamples, Lu(ii)=Lx(u(ii,:),dt,nodes.U(ii),'expoDecay'); Lv(ii)=Lx(v(ii,:),dt,nodes.U(ii),'expoDecay'); Lw(ii)=Lx(w(ii,:),dt,nodes.U(ii),'expoDecay'); end figure box on;hold on,grid on plot(nodes.Y,Lu,'ko','markerFaceColor','r') plot(nodes.Y,Lv,'ko','markerFaceColor','y') plot(nodes.Y,Lw,'ko','markerFaceColor','b') legend('Lu','Lv','Lw') ylim([0,220]) xlabel('y (m)') ylabel(' Integral length scales (m)') set(gcf,'color','w') % set(findall(gca,'-property','fontSize'),'fontsize',14) 

## Standard deviation of wind velocity fluctuations

figure box on;hold on,grid on plot(nodes.Y, std(u,0,2),'ko','markerFaceColor','r') plot(nodes.Y, std(v,0,2),'ko','markerFaceColor','y') plot(nodes.Y, std(w,0,2),'ko','markerFaceColor','b') ylim([0,3]) xlabel('y (m)') ylabel('\sigma_{i} (m/s), i = (u,v,w)') set(gcf,'color','w') legend('\sigma_{u}','\sigma_{v}','\sigma_{w}','location','SouthWest') % set(findall(gca,'-property','fontSize'),'fontsize',16) 

## Turbulence Intensity

Iu = std(u,0,2)./nodes.U*100; Iv = std(v,0,2)./nodes.U*100; Iw = std(w,0,2)./nodes.U*100; figure box on;hold on,grid on plot(nodes.Y,Iu,'ko','markerFaceColor','r') plot(nodes.Y,Iv,'ko','markerFaceColor','y') plot(nodes.Y,Iw,'ko','markerFaceColor','b') ylim([0,20]) xlabel('y (m)') ylabel('TI (%)') set(gcf,'color','w') legend('Iu','Iv','Iw','location','SouthWest') % set(findall(gca,'-property','fontSize'),'fontsize',16) 

## Comparison with Von Karman spectrum

block duration (each)

tBlock = 600; % number of blocks NBlock = round(tmax/tBlock); % number of data point per block Ncoh = round(N/NBlock); clear Su Sv Sw for ii=1:Nsamples, [Su(ii,:),f]=pwelch(u(ii,:),Ncoh,round(Ncoh/2),Ncoh,fs); [Sv(ii,:),f]=pwelch(v(ii,:),Ncoh,round(Ncoh/2),Ncoh,fs); [Sw(ii,:),f]=pwelch(w(ii,:),Ncoh,round(Ncoh/2),Ncoh,fs); end Su = mean(Su)'; % average the single-point spectra Sv = mean(Sv)'; % average the single-point spectra Sw = mean(Sw)'; % average the single-point spectra % empirical Von Karman spectra coef=[-5/6, -11/6]; Lu_target = 170; Lv_target = 100; Lw_target = 35; stdU_target = 2.6; stdV_target = 2.2; stdW_target = 1.56; U_target = 20; n = Lu_target./U_target.*f'; Su_target = U_target.^(-1).*4.*Lu_target.*stdU_target.^2.*(1+70.7.*n.^2).^(coef(1)); n = Lv_target./U_target.*f'; Sv_target= U_target.^(-1).*4.*Lv_target.*stdV_target.^2.*(1+70.7.*4.*n.^2).^(coef(2)).*(1+188.4.*4*n.^2); n = Lw_target./U_target.*f'; Sw_target= U_target.^(-1).*4.*Lw_target.*stdW_target.^2.*(1+70.7.*4.*n.^2).^(coef(2)).*(1+188.4.*4*n.^2); %plot the spectra figure box on;hold on plot(f,f.*Su_target','r') plot(f,f.*Su,'k') legend('Von Karman','computed') xlabel(' freq (Hz) '); set(gca,'Xscale','log') ylabel('f \cdot S_{u} ( m^2 \cdot s^{-2})') set(gcf,'color','w') xlim([1/600,fs/2]) figure box on;hold on plot(f,f.*Sv_target','r') plot(f,f.*Sv,'k') legend('Von Karman','computed') xlabel(' freq (Hz) '); set(gca,'Xscale','log') ylabel('f \cdot S_{v} ( m^2 \cdot s^{-2})') set(gcf,'color','w') xlim([1/600,fs/2]) figure box on;hold on plot(f,f.*Sw_target','r') plot(f,f.*Sw,'k') legend('Von Karman','computed') xlabel(' freq (Hz) '); set(gca,'Xscale','log') ylabel('f \cdot S_{w} ( m^2 \cdot s^{-2})') set(gcf,'color','w') xlim([1/600,fs/2]) 

## Co-coherence calcualtion

The co-coherence is preferred to the magnitude squared coehrence when it comes to wind flucutation block duration (each)

tBlock = 240; % number of blocks NBlock = round(tmax/tBlock); % number of data point per block Ncoh = round(N/NBlock); if mod(Ncoh,2), cocoh = zeros(Nsamples,Nsamples,round(Ncoh/2)); else cocoh = zeros(Nsamples,Nsamples,round(Ncoh/2)+1); end % computation of the cocoherence using the function coherence for ii=1:Nsamples, for jj=1:Nsamples, [cocoh(ii,jj,:),~,freq] = coherence(u(ii,:),u(jj,:),Ncoh,round(Ncoh/2),Ncoh,fs); end end % matrix distance along y dy = zeros(Nsamples,Nsamples); for kk= 1:Nsamples, for mm = 1:Nsamples, dy(kk,mm) = abs(nodes.Y(kk)-nodes.Y(mm)); end end % average of coherence over same cross-flow separations cocoh = reshape(cocoh,Nsamples*Nsamples,[]); [uniqueDist]=unique(round(dy(:)*100)/100); % average at +/- 1 cm meanCoCoh = zeros(numel(uniqueDist),size(cocoh,2)); for ii=1:numel(uniqueDist), ind = find(round(dy(:)*100)/100==uniqueDist(ii)); meanCoCoh(ii,:)=nanmean(cocoh(ind,:)); end 

## Co-coherence verification

% choice of 4 target distances distTarget = [10,20,40,60]; % wave number k k = 2*pi*freq./nodes.U(1); figure for ii=1:4, [~,indDist] = min(abs(uniqueDist-distTarget(ii))); subplot(2,2,ii) box on;hold on plot(k,meanCoCoh(indDist,:),'ko','markerfacecolor','y'); plot(k,exp(-7.*dy(1,indDist).*freq/nodes.U(1)),'k','linewidth',2); legend(['d = ',num2str(uniqueDist(indDist),2),' m'],'target'); xlim([0,0.12]) ylabel('co-coherence') if ii>2, xlabel('k (m^{-1})') end set(gcf,'color','w') end 

## Conclusions

A simple example of turbulent wind simulation is presented here. Because the simulation procedure belongs to Monte-Carlo schemes, the measured parameters are not always expected to perfectly fit with the target ones. At every new simulation, different time series are obtained. The user is free to implement new wind spectra or new coherence models.