help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Octave forge] [bim package] bim1a_advection*


From: JuanPi
Subject: [Octave forge] [bim package] bim1a_advection*
Date: Sat, 11 Jan 2020 00:56:27 +0100

Hi again,

I keep playing with the bim package.
I am failing to understand the behavior (or likely the correct usage)
of bim1a_advection_upwind.
I noticed that for finer meshes the advection effect decreases, and
converges to a solution with no advection at all.
How can doe sone create correct (almost mesh independent, for
resoanable meshes) advection-diffusion matrices?

** Example of the problem

Given this equation

Y_t - a Y_xx + b Y_x + Y = f(x)
Y(0, t) = Yo
Y_t(1,t) = -b Y_x(1,t)   # open boundary approx O(a^2)

with Y_<var> is the derivative of Y w.r.t. <var>

Compare the solutions generated with this code

a  = 1e-4;
b  = 100;
Yo = 1;

# Initial condition
Y0 = @(x) Yo * ones (size (x));

# Sources
c = 10;
f = @(x)Yo + (c - Yo) * exp(-(x-0.5).^2 / 2 / 0.05^2);

# Dynamic equations
function dY = dYdt(Y, t, b, A, msh, f)
 dY      = f - A * Y;
 dY(1)   = 0;
 dY(end) = - b * (Y(end) - Y(end-1)) / (msh(end) - msh(end-1));
endfunction
lsode_options ("integration method", "adams");

# Meshes
Nt     = 10;
t      = linspace (0, 10, Nt);
msh{1} = linspace (0, 1, 50).';
msh{2} = linspace (0, 1, 100).';

# Matrices & solution
Nmsh = numel(msh);
for i = 1:Nmsh
  Diff     = bim1a_laplacian (msh{i}, 1, a);
  Conv     = bim1a_advection_upwind (msh{i}, b * msh{i});
  Reaction = bim1a_reaction (msh{i}, 1, 1);
  Source   = bim1a_rhs (msh{i}, 1, f(msh{i}));

  A = Diff + Conv + Reaction;

  printf ('Solving %d\n', i); fflush (stdout)
  Y{i} = lsode (@(Y, t) dYdt(Y, t, b, A, msh{i}, Source), Y0(msh{i}), t);
endfor # over meshes

figure (1)
col = flipud (gray (Nt+3))(3:end,:); # colors for time
for i=1:Nmsh
  subplot (2, 1, i);
  h_ = plot (msh{i}, Y{i}.');

  for j=1:Nt
    set(h_(j), 'color', col(j,:));
  endfor # over colors

  ylabel (sprintf('Y%d',i))
endfor # over solutions
xlabel ('x')


Thanks

-- 
JuanPi Carbajal
https://goo.gl/ayiJzi

-----
“An article about computational result is advertising, not
scholarship. The actual scholarship is the full software environment,
code and data, that produced  the  result.”
- Buckheit and Donoho



reply via email to

[Prev in Thread] Current Thread [Next in Thread]