Summation function with two variables concentration profile

Hi,

Need some help plotting the attached equation For some reason my plot is not oming out as needed. For example, I want to be able to plot C vs X for this equation. For some reason my methodology doesn't com out right

D=1.79e-5; 
C1=0.79; 
C0=0; 
x=linspace(-60,60,1000);
n=linspace(0,.020,1000);
t=linspace(0,1,1000);
xs=linspace(0,120,100);

ax1=(4./pi()).*(((-1).^n)./(2.*n+1));
bx1=exp(-((D*((2*n+1).^2)*(pi().^2).*t)./(4.*l.^2)));

Unrecognized function or variable 'l'.

cx1=cos(((2.*n+1).*pi().*x)./(2.*l));
CC1=C0+(C1-C0).*(1-ax1.*bx1.*cx1);
figure(1)
plot(xs,CC1,'-k')

Answer 1

[Walter Roberson]

D=1.79e-5; 
C1=0.79; 
C0=0; 
l = 120/2;
x=linspace(-l,l,1000);
n=linspace(0,.020,1000).';
t=linspace(0,1,1000).';

ax1=(4./pi()).*(((-1).^n)./(2.*n+1));
bx1=exp(-((D*((2*n+1).^2)*(pi().^2).*t)./(4.*l.^2)));
cx1=cos(((2.*n+1).*pi().*x)./(2.*l));
CC1=C0+(C1-C0).*(1-ax1.*bx1.*cx1);
figure(1)
whos

  Name         Size                 Bytes  Class     Attributes

  C0           1x1                      8  double              
  C1           1x1                      8  double              
  CC1       1000x1000            16000000  double    complex   
  D            1x1                      8  double              
  ax1       1000x1                  16000  double    complex   
  bx1       1000x1                   8000  double              
  cx1       1000x1000             8000000  double              
  l            1x1                      8  double              
  n         1000x1                   8000  double              
  t         1000x1                   8000  double              
  x            1x1000                8000  double              

figure
surf(x, t, real(CC1), 'edgecolor', 'none')
title('real CC1')

figure
surf(x, t, imag(CC1), 'edgecolor', 'none')
title('imag CC1')

[mpz]

@Walter Roberson do you know how I output the 2-D plot for x vs CC1(real part) only without t? Also, if I wanted to use a different initial timestep of say 3000s, do I do it as below?

x=linspace(-l,l,3000);

n=linspace(0,.020,3000).';

t=linspace(0,1,3000).';

[Walter Roberson]

D=1.79e-5; 
C1=0.79; 
C0=0; 
l = 120/2;
x=linspace(-l,l,1000);
n=linspace(0,.020,1000).';
t=linspace(0,1,1000).';

ax1=(4./pi()).*(((-1).^n)./(2.*n+1));
bx1=exp(-((D*((2*n+1).^2)*(pi().^2).*t)./(4.*l.^2)));
cx1=cos(((2.*n+1).*pi().*x)./(2.*l));
CC1=C0+(C1-C0).*(1-ax1.*bx1.*cx1);
figure(1)
whos

  Name         Size                 Bytes  Class     Attributes

  C0           1x1                      8  double              
  C1           1x1                      8  double              
  CC1       1000x1000            16000000  double    complex   
  D            1x1                      8  double              
  ax1       1000x1                  16000  double    complex   
  bx1       1000x1                   8000  double              
  cx1       1000x1000             8000000  double              
  l            1x1                      8  double              
  n         1000x1                   8000  double              
  t         1000x1                   8000  double              
  x            1x1000                8000  double              

figure
plot(x, real(CC1))
title('real CC1')

figure
plot(x, imag(CC1))
title('imag CC1')

[Walter Roberson]

Remember, you are defining which is a surface defined over two independent variables. If you do not have time as an independent axis in your plot, then you should be creating one plot for each different timestep.

[mpz]

@Walter Roberson that is what I was getting but the 2D plot for that equations should look like below. Looking at it from the 3D plot it looks correct but converting to 2D to look like below doesn't. Do you have a fix for that?

[Walter Roberson]

n=linspace(0,.020,1000);

If that n is intended to be the same n as your equation, then note that n for that equation comes in in the form of . But that means that n should be integers 0, 1, 2, 3 ... And that makes a big difference, as then your is strictly integer powers which would give you real values rather than complex values.

[mpz]

yes that n is supposed to be intergers only. Sorry for the confusion. I don't know how to code that summation as it goes to inf as well.

[mpz]

@Walter Roberson the intend is to plot C which is the concentration for each step size. I have tried setting t=1000, t=3000 etc. but it doesn't come out as needed. The end result is something as below for each time step.

[Walter Roberson]

You can evaluate over a range of values of n and sum() afterwards.

But what precision do you need? Although the calculation does appear to converge (perhaps I should have checked more carefully) the third digit is still changing when you take a million terms -- the terms might be small but they add up in bulk.

[Walter Roberson]

Do you have a fix for that?

No, I do not have a fix for that. I do not think you are calculating the correct thing.

For a calculation like this, every independent variable that is being varied over should have its own axes. x should be a different axes than t which should be a different axes than n. You would typically calculate for a finite number of n and then sum() over that axes, giving you a result that has x and t and possibly other axes. Once you have the multi-dimensional output you would try to summarize it into plots.

Because your x and t are different continuous independent axes, your plot should not be trying to drop t on the output.