Academic Integrity: tutoring, explanations, and feedback — we don’t complete graded work or submit on a student’s behalf.

I am done with the code and even got the plots but I can\'t find the right funct

ID: 3551238 • Letter: I

Question

I am done with the code and even got the plots but I can't find the right function to search for the peak concentrations for stability class B and D. I need the magnitude and the location as stated in part 2 of the problem.


I am attaching my code here as well. It's on Matlab.



% input variables

Q=500*10^6; %source emission rate (ug/s)

H=120; %effective stack height (m)


% define spatial domain

xlength=10000; % domain length (m) --> 10km

ylength=4000; % domain width (m) --> 2 km on either side of the plume centerline

dx=100; % x node spacing (m)

dy=100; % y node spacing (m)

x=0:dx:xlength; %(m)

x(1)=1;

y=-ylength/2:dy:ylength/2; %(m)

% y' or x' for transpose


% calculate windspeed at z=H for stability classes B, D

% ground-level windspeeds (at z=10 m)

GLwindspeed=[2 2]; %(m/s)



% uH/ua=(H/za)^p, so uH=ua*(H/za)^p. ua=anemometer height=GLwindspeed

% exponent p varies with stability class


p=[0.15 0.25]; % from Table 7.6

za=10; % standard anemometer height, (m)

% multiply the P values with 0.6 for smooth terrain


for i=1:2 % # of stability classes to be considered

uH(i)=GLwindspeed(i)*(H/za)^p(i); %wind speed at effective stack height (m)

end

%for two stability classes we will only have 2 values for p and i will go from 1:2



% coefficients a, b, d, f for calculation of sigma y and z (Table 7.8)

coeff_a=[156 68]; % coefficient a doesn't depend on the downwind distance x

coeff_c=[106.6 33.2; 108.2 44.5]; % different coefficients for x<=1km and x>1km

coeff_d=[1.149 0.725; 1.098 0.516]; % different coefficients for x<=1km and x>1km

coeff_f=[3.3 -1.7; 2 -13]; % different coefficients for x<=1km and x>1km


%sigma y = ax^0.894, sigma z= cx^d+f


%calculate the pollutant concentration C for entire domain (and the two stability classes)


for st=1:2 % st = # of stability class (B and D)

  

for i=1:length(x)

  

sigmaY=coeff_a(st)*(x(i)/1000)^0.894; % sigma Y doesn't depend on the downwind distance x

  

% sigma Z depends on the downwind distance x:

if x(i)<=1000 %(m)

sigmaZ=coeff_c(1,st)*(x(i)/1000)^coeff_d(1,st)+coeff_f(1,st);

else

sigmaZ=coeff_c(2,st)*(x(i)/1000)^coeff_d(2,st)+coeff_f(2,st);

end

  

for j=1:length(y)

  

C(i,j,st)=Q/(pi*sigmaY*sigmaZ*uH(st))*exp(-y(j)^2/(2*sigmaY^2))*2*exp(-H^2/(2*sigmaZ^2));

  

end

end

end




please help me figure out how to search for the max concentrations.

In this project, you will use the Point Source Gaussian Plume Model to approximate concentrations of a pollutant downwind of a power plant. The Point Source Gaussian Plume Model is given by where x is distance along the plume centerline (m), y is the horizontal distance from the plume centerline (m), z is the vertical distance from ground surface at the emission site (m), Q is the source emission rate (mug/s), H is the effective stack height (m), uH is the average wind velocity at z = H, sigmay and sigmaz are the dispersion coefficients in the y and z directions (m). The average emission rate of SO2 from the power plant is 500 g/s, from a stack of height 100 m (assume H = 120 m). The topography in the region is relatively flat so we can assume a ground surface at z=0 everywhere. The ground-level windspeed is 2m/s (for both stability classes B and D). Use the MATLAB code developed in class to calculate the ground-level S02 concentration for 10 km downwind of the stack and 2 km on either side of the plume centerline (a 40 km2 area), for stability classes B and D. Generate a contour plot of ground-level concentration for each scenario. Find the location and magnitude of peak downwind concentrations for each stability class). Explain your results: how do the location and magnitude of the peak concentration vary with stability class and why? Hint use Matlab functions max"and m[ij]=find0".

Explanation / Answer

% input variables

Q=500*10^6; %source emission rate (ug/s)

H=120; %effective stack height (m)


% define spatial domain

xlength=10000; % domain length (m) --> 10km

ylength=4000;   % domain width (m) --> 2 km on either side of the plume centerline

dx=100;         % x node spacing (m)

dy=100;         % y node spacing (m)

x=0:dx:xlength; %(m)

x(1)=1;

y=-ylength/2:dy:ylength/2; %(m)

% y' or x' for transpose


% calculate windspeed at z=H for stability classes B, D

% ground-level windspeeds (at z=10 m)

GLwindspeed=[2 2]; %(m/s)



% uH/ua=(H/za)^p, so uH=ua*(H/za)^p. ua=anemometer height=GLwindspeed

% exponent p varies with stability class


p=[0.15 0.25]; % from Table 7.6

za=10; % standard anemometer height, (m)

% multiply the P values with 0.6 for smooth terrain


for i=1:2 % # of stability classes to be considered

    uH(i)=GLwindspeed(i)*(H/za)^p(i); %wind speed at effective stack height (m)

end

%for two stability classes we will only have 2 values for p and i will go from 1:2



% coefficients a, b, d, f for calculation of sigma y and z (Table 7.8)

coeff_a=[156 68]; % coefficient a doesn't depend on the downwind distance x

coeff_c=[106.6 33.2; 108.2 44.5]; % different coefficients for x<=1km and x>1km

coeff_d=[1.149 0.725; 1.098 0.516]; % different coefficients for x<=1km and x>1km

coeff_f=[3.3 -1.7; 2 -13]; % different coefficients for x<=1km and x>1km


%sigma y = ax^0.894, sigma z= cx^d+f


%calculate the pollutant concentration C for entire domain (and the two stability classes)


for st=1:2 % st = # of stability class (B and D)

   

    for i=1:length(x)

       

        sigmaY=coeff_a(st)*(x(i)/1000)^0.894; % sigma Y doesn't depend on the downwind distance x

       

        % sigma Z depends on the downwind distance x:

      

         if x(i)<=1000 %(m)

            sigmaZ=coeff_c(1,st)*(x(i)/1000)^coeff_d(1,st)+coeff_f(1,st);

        else

            sigmaZ=coeff_c(2,st)*(x(i)/1000)^coeff_d(2,st)+coeff_f(2,st);

        end

       

        for j=1:length(y)

           

            C(i,j,st)=Q/(pi*sigmaY*sigmaZ*uH(st))*exp(-y(j)^2/(2*sigmaY^2))*2*exp(-H^2/(2*sigmaZ^2));

           

        end

    end

end
% for stability class B
max_c=max(max(C(:,:,1)))% apply max two times since C is 2D array
[i,j]=find(C(:,:,1)==max_c)% i,j have the value where the peaks occur
% for stability class D
max_c=max(max(C(:,:,2)))
[i,j]=find(C(:,:,2)==max_c)