Added additional boundary conditions

Source Link

edited May 22, 2018 at 16:21

43
4

FiniteElement and NDSolve::bdord error won't work

I need to solve a system of DE's in 2 dimensions over time, The problem is that mathematica complains when I have a boundary condition of the same order as the differential order of the DE's:

NDSolve::bdord: Boundary condition (u^(0,1,0))[x,0,t] should have derivatives of order lower than the differential order of the partial differential equation.

I'm curious if there are any methods regarding NDSolve that allow, so I tried using the finiteelement method since it worked for thisother people who had similar problems, or does this need to be done another way?but it just returns itself as if nothing ever happened

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t], 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t], 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

(D[u[x, y, t], x] /. x -> 0) == (D[u[x, y, t], x] /. x -> 12000000),
(D[v[x, y, t], x] /. x -> 0) == (D[v[x, y, t], x] /. x -> 12000000),
D[u[0, y, t], y] == D[u[12000000, y, t], y],
D[v[0, y, t], y] == D[v[12000000, y, t], y],
D[u[0, y, t], t] == D[u[12000000, y, t], t],
D[v[0, y, t], t] == D[v[12000000, y, t], t],

(D[T[x, y, t], x] /. x -> 0) == (D[T[x, y, t], x] /. x -> 12000000),
(D[ρ[x, y, t], x] /. x -> 0) == (D[ρ[x, y, t], x] /. 
x -> 12000000),
D[T[0, y, t], y] == D[T[12000000, y, t], y],
D[ρ[0, y, t], y] == D[ρ[12000000, y, t], y],
D[T[0, y, t], t] == D[T[12000000, y, t], t],
D[ρ[0, y, t], t] == D[ρ[12000000, y, t], t],

(D[e[x, y, t], x] /. x -> 0) == (D[e[x, y, t], x] /. x -> 12000000),
D[e[0, y, t], y] == D[e[12000000, y, t], y],
D[e[0, y, t], t] == D[e[12000000, y, t], t],


e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100},
Method -> "FiniteElement" ]

NDSolve::bdord error

I need to solve a system of DE's in 2 dimensions over time, The problem is that mathematica complains when I have a boundary condition of the same order as the differential order of the DE's:

NDSolve::bdord: Boundary condition (u^(0,1,0))[x,0,t] should have derivatives of order lower than the differential order of the partial differential equation.

I'm curious if there are any methods regarding NDSolve that allow for this, or does this need to be done another way?

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t], 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t], 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

FiniteElement and NDSolve won't work

I need to solve a system of DE's in 2 dimensions over time, The problem is that mathematica complains when I have a boundary condition of the same order as the differential order of the DE's, so I tried using the finiteelement method since it worked for other people who had similar problems, but it just returns itself as if nothing ever happened

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t], 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t], 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

(D[u[x, y, t], x] /. x -> 0) == (D[u[x, y, t], x] /. x -> 12000000),
(D[v[x, y, t], x] /. x -> 0) == (D[v[x, y, t], x] /. x -> 12000000),
D[u[0, y, t], y] == D[u[12000000, y, t], y],
D[v[0, y, t], y] == D[v[12000000, y, t], y],
D[u[0, y, t], t] == D[u[12000000, y, t], t],
D[v[0, y, t], t] == D[v[12000000, y, t], t],

(D[T[x, y, t], x] /. x -> 0) == (D[T[x, y, t], x] /. x -> 12000000),
(D[ρ[x, y, t], x] /. x -> 0) == (D[ρ[x, y, t], x] /. 
x -> 12000000),
D[T[0, y, t], y] == D[T[12000000, y, t], y],
D[ρ[0, y, t], y] == D[ρ[12000000, y, t], y],
D[T[0, y, t], t] == D[T[12000000, y, t], t],
D[ρ[0, y, t], t] == D[ρ[12000000, y, t], t],

(D[e[x, y, t], x] /. x -> 0) == (D[e[x, y, t], x] /. x -> 12000000),
D[e[0, y, t], y] == D[e[12000000, y, t], y],
D[e[0, y, t], t] == D[e[12000000, y, t], t],


e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100},
Method -> "FiniteElement" ]

Changed the problem, since I cannot use second order derivatives, the problem just becomes one regarding the NDSolve::bdord

Source Link

edited May 22, 2018 at 13:45

Hidde Rinsema

43
4

Stiff system error in NDSolve::bdord error

I haveneed to solve a stiff system errorof DE's in 2 dimensions over time, any solution The problem is that mathematica complains when I have tried tries to evaluate ita boundary condition of the same order as the differential order of the DE's:

NDSolve::bdord: Boundary condition (u^(0,1,0))[x,0,t] should have derivatives of order lower than the differential order of the partial differential equation.

I'm curious if there are any methods regarding NDSolve that allow for this, but then stops(don't know why).or does this need to be done another way?

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
η = 10^-6;
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t] - η (D[u[x, y, t], y] + D[v[x, y, t], x]), 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t] - η (D[v[x, y, t], x] + D[u[x, y, t], x]), 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

The error is as follows:

NDSolve::ndsz: At t == 0.0005560756736763107`, step size is effectively   zero; singularity or stiff system suspected.

I also get this error:

NDSolve::eerr: Warning: scaled local spatial error estimate of 332.45447642223894` at t = 0.0005560756736763107` in the direction of independent variable x is much greater than the prescribed error tolerance. Grid spacing with 13 points may be too large to achieve the desired accuracy or precision. A singularity may have formed or a smaller grid spacing can be specified using the MaxStepSize or MinPoints method options.

I do not know if the first error causes the second or visa-versa, my area is rather large, but I do not know how to fix the gridspace problem, as I said, when I try to solve it with possible solutions(increasing maximum stepsize for example) the thing stops for no given reason, perhaps this has something to do with ram or maximum running time. As you can see, I'm not that experienced with solving stuff in higher than 1 dimension, I hope this is a resolvable problem.

EDIT: When I add a MaxStepFraction of 1/100 the code runs, but the cell keeps running and doesn't produce results. But unlike my comments earlier it keeps running and doesn't stop randomly

Stiff system error in NDSolve

I have a stiff system error, any solution I have tried tries to evaluate it, but then stops(don't know why).

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
η = 10^-6;
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t] - η (D[u[x, y, t], y] + D[v[x, y, t], x]), 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t] - η (D[v[x, y, t], x] + D[u[x, y, t], x]), 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

The error is as follows:

NDSolve::ndsz: At t == 0.0005560756736763107`, step size is effectively   zero; singularity or stiff system suspected.

I also get this error:

NDSolve::eerr: Warning: scaled local spatial error estimate of 332.45447642223894` at t = 0.0005560756736763107` in the direction of independent variable x is much greater than the prescribed error tolerance. Grid spacing with 13 points may be too large to achieve the desired accuracy or precision. A singularity may have formed or a smaller grid spacing can be specified using the MaxStepSize or MinPoints method options.

I do not know if the first error causes the second or visa-versa, my area is rather large, but I do not know how to fix the gridspace problem, as I said, when I try to solve it with possible solutions(increasing maximum stepsize for example) the thing stops for no given reason, perhaps this has something to do with ram or maximum running time. As you can see, I'm not that experienced with solving stuff in higher than 1 dimension, I hope this is a resolvable problem.

EDIT: When I add a MaxStepFraction of 1/100 the code runs, but the cell keeps running and doesn't produce results. But unlike my comments earlier it keeps running and doesn't stop randomly

NDSolve::bdord error

I need to solve a system of DE's in 2 dimensions over time, The problem is that mathematica complains when I have a boundary condition of the same order as the differential order of the DE's:

NDSolve::bdord: Boundary condition (u^(0,1,0))[x,0,t] should have derivatives of order lower than the differential order of the partial differential equation.

I'm curious if there are any methods regarding NDSolve that allow for this, or does this need to be done another way?

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t], 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t], 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

Mistake in that t wasn't added for e

Source Link

edited May 22, 2018 at 11:20

Hidde Rinsema

43
4

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu );
 kb T[x, y, t];
η = 10^-6;
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t] - η (D[u[x, y, t], y] + D[v[x, y, t], x]), 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t] - η (D[v[x, y, t], x] + D[u[x, y, t], x]), 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu );
kb T[x, y, t];
η = 10^-6;
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t] - η (D[u[x, y, t], y] + D[v[x, y, t], x]), 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t] - η (D[v[x, y, t], x] + D[u[x, y, t], x]), 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

P[x_, y_, t_] = e[x, y, t]/(γ - 1) ; 
e[x_, y_, t_] = (γ - 1) ρ[x, y, t]/(μ mu ) kb T[x, y, t];
η = 10^-6;
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;
sol1 = NDSolve[{
D[ρ[x, y, t]*u[x, y, t], 
 t] == -D[ρ[x, y, t]*u[x, y, t]*u[x, y, t] + P[x, y, t], x] -
  D[ρ[x, y, t]*u[x, y, t]*
    v[x, y, t] - η (D[u[x, y, t], y] + D[v[x, y, t], x]), 
  y],
D[ρ[x, y, t]*v[x, y, t], 
 t] == -D[ρ[x, y, t]*v[x, y, t]*
     u[x, y, t] - η (D[v[x, y, t], x] + D[u[x, y, t], x]), 
   y] - D[ρ[x, y, t]*v[x, y, t]*v[x, y, t] + 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]),

v[0, y, t] == v[12*10^6, y, t],
u[0, y, t] == u[12*10^6, y, t],
T[0, y, t] == T[12*10^6, y, t],
ρ[0, y, t] == ρ[12*10^6, y, t],

e[x, 0, t] == 3.83767261162,
v[x, 4000000, t] == 0,
v[x, 0, t] == 0,
(D[u[x, y, t], y] /. y -> 0) == 0,
(D[u[x, y, t], y] /. y -> 4000000) == 0,

v[x, y, 0] == 0,
u[x, y, 0] == 0,
T[x, y, 0] == 5770 + 0.00835414960707927 y,
ρ[x, y, 0] == 
1.42*10^-7*1.408*10^3 + 7.3561137493644*10^-10 y
},
{u, v, T, ρ}, {x, 0, 12000000}, {y, 0, 4000000}, {t, 0, 100}]

deleted 260 characters in body

Source Link

edited May 22, 2018 at 9:37

Hidde Rinsema

43
4

Loading

Added information regarding the addition of MaxStepFraction

Source Link

edited May 21, 2018 at 14:21

Hidde Rinsema

43
4

Loading

added 1 character in body

Source Link

edited May 21, 2018 at 13:57

Hidde Rinsema

43
4

Loading

deleted 147 characters in body

Source Link

edited May 21, 2018 at 12:24

Henrik Schumacher

109.5k
7
186
323

Loading

Source Link

asked May 21, 2018 at 12:15

Hidde Rinsema

43
4

Loading

Stack Exchange Network

Return to Question

FiniteElement and NDSolve::bdord error won't work

NDSolve::bdord error

FiniteElement and NDSolve won't work

Stiff system error in NDSolve::bdord error

Stiff system error in NDSolve

NDSolve::bdord error