%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Exercise set 7: Properties of turbulence %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Instructions: Complete the following matlab script to accomplish the % tasks required by the exercise set. Much of the script is done for you, % you will need to complete lines marked with 'ReplaceMe' %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clc; % Clear the command window clear; % Clear all stored variables close all; % Close all open figures ReplaceMe = NaN; % DO NOT CHANGE THIS VALUE % Import the data %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% u = load('u.out'); % Load u-velocity data w = load('w.out'); % Load w-velocity data q = load('q.out'); % Load specific humidity data % ACF function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % please look the end of the file %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % a. Plot u,w, and q as a function of time. % % Mark the mean of u w and q on each plot. % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Step 1. Create a time vector, with the same size as the data, % corresponding to a sampling frequency of 20 Hz. fs = ReplaceMe; % Sampling frequency [Hz] N = ReplaceMe; % Number of samples t = ReplaceMe; % Time vector % Step 2. Determine the mean values of the data u_avg = ReplaceMe; % Find the mean value of u; w_avg = ReplaceMe; % Find the mean value of w; q_avg = ReplaceMe; % Find the mean value of q; % Print the result to the command window fprintf('Mean of u = %3.2f [m/s] \n', u_avg) fprintf('Mean of w = %3.2e [m/s] \n', w_avg) fprintf('Mean of q = %3.2f [g/m^3] \n',q_avg) fprintf('Total sampling time = %3.2e [s] \n\n',t(N)) % Step 3. Plot the raw data in terms of time, % Create a figure that is 16x16cm with a white background % plots will appear here. figure('color','w','unit','centimeters','position',[1 4 16 16]) ax{1} = subplot(311); % Create a subplot, of a 3x1 set of subplots, % this will be the first, (hence, 311). Save a % reference to the axis to ax{1} ax{2} = subplot(312); % Create a subplot, of a 3x1 set of subplots, % this will be the second, (hence, 312). Save a % reference to the axis to ax{2} ax{3} = subplot(313); % Create a subplot, of a 3x1 set of subplots, % this will be the third, (hence, 313). Save a % reference to the axis to ax{3} hold(ax{1},'on') % Allow multiple plots on the same axes hold(ax{2},'on') % Allow multiple plots on the same axes hold(ax{3},'on') % Allow multiple plots on the same axes box(ax{1},'on') % Display the axes as a box box(ax{2},'on') % Display the axes as a box box(ax{3},'on') % Display the axes as a box ylabel(ax{1},'u [m/s]') % Label the y-axis ylabel(ax{2},'w [m/s]') % Label the y-axis ylabel(ax{3},'q [g/m^3]') % Label the y-axis xlabel(ax{1},'t [s]') % Label the x-axis xlabel(ax{2},'t [s]') % Label the x-axis xlabel(ax{3},'t [s]') % Label the x-axis [ymin, ymax] = bounds([u w q]); % Find the upper and lower bounds of u [xmin, xmax] = bounds(t); % Find the upper and lower bounds of t axis(ax{1},[xmin xmax ymin(1) ymax(1)]); % Adjust the size of the plot axes; axis(ax{2},[xmin xmax ymin(2) ymax(2)]); % Adjust the size of the plot axes; axis(ax{3},[xmin xmax ymin(3) ymax(3)]); % Adjust the size of the plot axes; plot(ax{1},ReplaceMe,ReplaceMe,'b-'); % Plot u as a function of time % using a solid blue line 'b-' plot(ax{2},ReplaceMe,ReplaceMe,'b-'); % Plot w as a function of time % using a solid blue line 'b-' plot(ax{3},ReplaceMe,ReplaceMe,'b-'); % Plot q as a function of time % using a solid blue line 'b-' plot(ax{1},ReplaceMe,ReplaceMe,'r--'); % Plot u_avg on the same % axes as the u-time series plot(ax{2},ReplaceMe,ReplaceMe,'r--'); % Plot w_avg on the same % axes as the w-time series plot(ax{3},ReplaceMe,ReplaceMe,'r--'); % Plot q_avg on the same % axes as the q-time series legend(ax{1},'u(t)','u_{avg}') % Label the time series and mean values legend(ax{2},'w(t)','w_{avg}') % Label the time series and mean values legend(ax{3},'q(t)','q_{avg}') % Label the time series and mean values %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % b. Verify that the average of u', w', and q' = 0 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Extract the fluctuating component of each time series, based on % the principle u(t) = u_avg + u', where u' is the fluctuating component u_fluct = ReplaceMe; w_fluct = ReplaceMe; q_fluct = ReplaceMe; % Take the mean of the fluctuating component of each time series u_fluct_avg = ReplaceMe; w_fluct_avg = ReplaceMe; q_fluct_avg = ReplaceMe; % Print the result to the command window fprintf('Mean of u'' = %3.2e [m/s] \n', u_fluct_avg) fprintf('Mean of w'' = %3.2e [m/s] \n', w_fluct_avg) fprintf('Mean of q'' = %3.2e [g/m^3] \n\n',q_fluct_avg) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % c. Plot u', w', and q' and calculate the standard deviation % % of each times series % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Step 1. Plot the time series of u', w', and q' % Create a figure that is 16x16cm with a white background % plots will appear here. figure('color','w','unit','centimeters','position',[1 4 16 16]) ax{1} = subplot(311); % Create a subplot, of a 3x1 set of subplots, % this will be the first, (hence, 311). Save a % reference to the axis to ax{1} ax{2} = subplot(312); % Create a subplot, of a 3x1 set of subplots, % this will be the second, (hence, 312). Save a % reference to the axis to ax{2} ax{3} = subplot(313); % Create a subplot, of a 3x1 set of subplots, % this will be the third, (hence, 313). Save a % reference to the axis to ax{3} box(ax{1},'on') % Display the axes as a box box(ax{2},'on') % Display the axes as a box box(ax{3},'on') % Display the axes as a box [ymin, ymax] = bounds([u_fluct w_fluct q_fluct]); % Find the upper and lower bounds of u [xmin, xmax] = bounds(t); % Find the upper and lower bounds of t axis(ax{1},[xmin xmax ymin(1) ymax(1)]); % Adjust the size of the plot axes; axis(ax{2},[xmin xmax ymin(2) ymax(2)]); % Adjust the size of the plot axes; axis(ax{3},[xmin xmax ymin(3) ymax(3)]); % Adjust the size of the plot axes; plot(ax{1},ReplaceMe,ReplaceMe,'b-'); % Plot u' as a function of time % using a solid blue line 'b-' plot(ax{2},ReplaceMe,ReplaceMe,'b-'); % Plot w' as a function of time % using a solid blue line 'b-' plot(ax{3},ReplaceMe,ReplaceMe,'b-'); % Plot q' as a function of time % using a solid blue line 'b-' ylabel(ax{1},'u'' [m/s]') % Label the y-axis ylabel(ax{2},'w'' [m/s]') % Label the y-axis ylabel(ax{3},'q'' [g/m^3]') % Label the y-axis xlabel(ax{1},'t [s]') % Label the x-axis xlabel(ax{2},'t [s]') % Label the x-axis xlabel(ax{3},'t [s]') % Label the x-axis % Step 2. Calculate the standard deviation u_std = ReplaceMe; % Calculate the standard deviation of u w_std = ReplaceMe; % Calculate the standard deviation of w q_std = ReplaceMe; % Calculate the standard deviation of q % Print results to the command window fprintf('Standard deviation of u = %3.2e [m/s] \n', u_std) fprintf('Standard deviation of w = %3.2e [m/s] \n', w_std) fprintf('Standard deviation of q = %3.2e [m/s] \n\n', q_std) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % d. Determine the turbulent fluxes of momentum u'w'_avg % % and vapor, w'q'_avg % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% uw_avg = ReplaceMe; % Determine the turbulent flux of momentum wq_avg = ReplaceMe; % Determine the turbulent flux of vapor % Print the result to the command window fprintf('Mean of u''w'' = %3.2e [m^2/s^2] \n', uw_avg) fprintf('Mean of w''q'' = %3.2e [g m^-2 s^-1] \n\n', wq_avg) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % e. based on the result of (d.), what time of day were these % % measurements taken? % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % f. Knowing the vertical gradient of the mean water vapor % % concentration at the measurement height (dqdz = -0.5 g m^-4)% % obtain the eddy diffusion coefficient for water vapor. % % Compare the value you find to the molecular diffusion % % coefficient Dm = 2.42e-5 m^2/s @ T=20C % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% dqdz= -0.5; % Vertical gradient of the mean water vapor concentration at % the measurement height Dt= ReplaceMe; fprintf('Eddy diffusion coefficient for water vapor: Dt = %3.2e [m^2/s] \n\n', Dt) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % g. Compute and plot the aurtocorrelation function of u', % % normalized so the maximum correlation is 1.0. you can use % % MATLAB's built-in autocorr function for this % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% u_acf = ReplaceMe; % Compute the complete autocorrelation function of u(t) % Step 1. Plot the time series of u', w', and q' % Create a figure that is 16x16cm with a white background % plots will appear here. figure('color','w','unit','centimeters','position',[1 4 16 16]) ax = axes(gcf); % Create an axes reference box(ax,'on') % Display the axes as a box [ymin, ymax] = bounds(u_acf); % Find the upper and lower bounds of u [xmin, xmax] = bounds(t); % Find the upper and lower bounds of t axis(ax,[xmin xmax ymin ymax]); % Adjust the size of the plot axes; plot(ax,ReplaceMe,ReplaceMe,'b-'); % Plot autocorrelation function as a function of % time-lags using a solid blue line 'b-' ylabel(ax,'\rho') % Label the y-axis xlabel(ax,'\tau [s]') % Label the x-axis %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % h. Calculate the integral time scale by integrating the % % autocorrelation function between tau = 0 and when the acf % % first becomes 0 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% tau_f = ReplaceMe; % Time when the autocorrelation function first becomes 0; Ti = ReplaceMe; % Integrate over dtau from 0