0
$\begingroup$

I use the following code to try and model 2D convection in time,

P = ρ[x, y, t]/(μ mu) kb T[x, y, t]  ;
e = P/(γ - 1);
cp = 5/2 kb/(μ mu);
Rgas = 8.3144598;
cv = 5/2 kb/(μ mu) - Rgas;
γ = cp/cv;
g = 28.02*9.81;
μ = 0.6163328197226503`;
mu = 1.66053904*10^-27;
kb = 1.38064852*10^-23;
NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
t] == -D[ρ[x, y, t]*(u[x, y, t])^2 + P, x],
D[ρ[x, y, t]*v[x, y, t], 
t] == -D[ρ[x, y, t]*(v[x, y, t])^2 + P, y] + 
g ρ[x, y, t],
D[ρ[x, y, t], t] == -D[ρ[x, y, t]*u[x, y, t], x] - 
D[ρ[x, y, t]*v[x, y, t], y],
D[e[x, y, t], t] == -D[u[x, y, t]*e[x, y, t], x] - 
D[v[x, y, t]*e[x, y, t], y] - 
P[x, y, t]*(D[u[x, y, t], x] - D[v[x, y, t], y]),

v[0, y, t] == v[12000000, y, t], u[0, y, t] == u[12000000, y, t], 
T[0, y, t] == 
T[12000000, y, t], ρ[0, y, t] == ρ[12000000, y, t],
D[u[x, y, t], y] /. y -> 4000000 == 0, 
D[u[x, y, t], y] /. y -> 0 == 0, v[x, 4000000, t] == 0, 
v[x, 0, t] == 0,
D[T[x, y, 0], y] == 0.41 g ρ [x, y, 0] T[x, y, 0]/P[x, y, 0], 
D[P[x, y, 0], y] == g ρ[x, y, 0], 
T[x, 0, 0] == 5770, ρ[x, 0, 0] == 1.42*10^-7*1.408*10^3
}, {u, v, 
T, ρ}, {{x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}}]

When evaluating this I get the following error

    NDSolve::dsvar: 0 cannot be used as a variable.   

Is that because I use u[x,y,0] and NDSolve sees the zero as a variable or are there other problems?

$\endgroup$
11
  • $\begingroup$ Your last two initial conditions, T[x, 0, 0] == 5770, \[Rho][x, 0, 0] == 1.42*10^-7*1.408*10^3 specify two variables to be zero, not just one. $\endgroup$
    – SPPearce
    Commented May 17, 2018 at 14:58
  • $\begingroup$ Oh, and your two initial conditions on $u_y$ need a parenthesis around the replacement rule. $\endgroup$
    – SPPearce
    Commented May 17, 2018 at 15:02
  • $\begingroup$ And you define P with no arguments, and then later try and give it three arguments. $\endgroup$
    – SPPearce
    Commented May 17, 2018 at 15:02
  • 1
    $\begingroup$ Oh, and finally the arguments at the end should not be in a list, this is what is causing the specific error message. $\endgroup$
    – SPPearce
    Commented May 17, 2018 at 15:25
  • 1
    $\begingroup$ @KraZug Would you like to sum up your comments in an answer? Or was it a simple mistake? $\endgroup$
    – Kuba
    Commented May 18, 2018 at 12:16

1 Answer 1

1
$\begingroup$

After fixing up the obvious typos, you get to:

cp = 5/2 kb/(μ mu);
Rgas = 8.3144598;
cv = 5/2 kb/(μ mu) - Rgas;
γ = cp/cv;
g = 28.02*9.81;
μ = 0.6163328197226503`;
mu = 1.66053904*10^-27;
kb = 1.38064852*10^-23;
P[x_, y_, t_] = ρ[x, y, t]/(μ mu) kb T[x, y, t];
e[x_, y_, t_] = P[x, y, t]/(γ - 1);

eqns = {D[ρ[x, y, t]*u[x, y, t], t] == -D[ρ[x, y, t]*(u[x, y, t])^2 + P[x, y, t], x],  
    D[ρ[x, y, t]*v[x, y, t], t] == -D[ρ[x, y, t]*(v[x, y, t])^2 + P[x, y, t], y] +  g ρ[x, y, t], 
    D[ρ[x, y, t], t] == -D[ρ[x, y, t]*u[x, y, t], x] - D[ρ[x, y, t]*v[x, y, t], y], 
    D[e[x, y, t], t] == -D[u[x, y, t]*e[x, y, t], x] - D[v[x, y, t]*e[x, y, t], y] - P[x, y, t]*(D[u[x, y, t], x] - D[v[x, y, t], y])}

bcs = {v[0, y, t] == v[12000000, y, t], u[0, y, t] == u[12000000, y, t], 
       T[0, y, t] == T[12000000, y, t], ρ[0, y, t] == ρ[12000000, y, t], 
       (D[u[x, y, t], y] /. y -> 4000000) == 0, (D[u[x, y, t], y] /. y -> 0) == 0, 
       v[x, 4000000, t] == 0, v[x, 0, t] == 0}

ics= {D[T[x, y, 0], y] == 0.41 g ρ[x, y, 0] T[x, y, 0]/P[x, y, 0], 
      D[P[x, y, 0], y] == g ρ[x, y, 0], T[x, y, 0] == 5770, 
      ρ[x, y, 0] == 1.42*10^-7*1.408*10^3}

NDSolve[Join[{eqns,bcs,ics}, {u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

However, this then gives an error to say that you can't include $u_y$ in the boundary conditions as the equations don't include a higher order derivative in $y$ (i.e. $u_{yy}$). If you remove these derivatives, it gives a different error. I suspect that you need to tell Mathematica explicitly which method you want it to use to discretize your equations, and it would be very sensible to nondimensionalise properly (i.e. don't have your $x$ and $y$ ranges go up to the millions).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.