Exploring Reyleigh-Taylor Instability with Variable Density and Viscosity

[MrOne]

Homework Statement

import time
import numpy as np
from Output import *
from PoissonSolver import *
tmps1=time.clock()

# Physical parameters
Re = 10.

# Numerical parameters
CFL = 0.9
dt = 0.001/4
t = 0.0
tmax = 1.5

# Make mesh
from MeshGeneration import *
mesh = Mesh(2.0,1.0,50,25)

# Output parameters
dtAff = 0.01 # s

# Boundary condition
# Cavity
#BC = {"type":{"Left": "Wall", "Right": "Wall","Bottom": "Wall", "Top": "MovingWall"},
# "u":{"Top": 1.0},
# "v":{}}
# Forcing Vortex
# BC = {"type":{"Left": "MovingWall", "Right": "MovingWall","Bottom": "MovingWall", "Top": "MovingWall"},
# "u":{"Top": 1.0, "Bottom": -1.0},
# "v":{"Left": 1.0, "Right": -1.0}}
# Channel
BC = {"type":{"Left": "Inflow", "Right": "Outflow","Bottom": "Wall", "Top": "Wall"},
"u":{"Left": 1.0},
"v":{"Left": 0.0},
"p":{"Right": 0.0}}
#BC = {"type":{"Left": "Wall", "Right": "Wall","Bottom": "Outflow", "Top": "Inflow"},
# "u":{"Top": 0.0},
# "v":{"Top": -1.0}}

# Initial condition
from InitNS import *
u, v, p = InitRest(mesh)

if BC["type"]["Left"] == "Inflow":
u[0,:] = BC["u"]["Left"]
if BC["type"]["Top"] == "Inflow":
v[:,mesh.Ny] = BC["v"]["Top"]

# Define fluxes
from NS_fluxes import GRAD, STRESS, TRANS, DIV

# Create Matrix for Poisson Solvre
A, bp = CreatePoissonMatrix(mesh, BC)

savefigNS(u,v,p,t,mesh)

while t<tmax:
ustar, vstar = np.copy(u), np.copy(v)
utrans, vtrans = TRANS(u,v,mesh,BC)
ustress, vstress = STRESS(u,v,mesh,BC)

ustar[1:-1,:] = u[1:-1,:] + dt*(utrans+1.0/Re*ustress)
vstar[:,1:-1] = v[:,1:-1] + dt*(vtrans+1.0/Re*vstress)
b = DIV(ustar,vstar,mesh)
p = PoissonSolver(A,b,bp,mesh)
#p = PoissonSolver(b,BC,mesh)
uprime, vprime = GRAD(p,mesh)
u, v = ustar-uprime, vstar-vprime

t = t + dt
# Plot
if (abs(t % dtAff) < dt):
savefigNS(u,v,p,t,mesh)
tmps2=time.clock()
print "Temps d'execution = %d\n" %(tmps2-tmps1)

Homework Equations

Hello everyone, i need your help to modify this program to include an scalar phi corresponding to a liquid fraction. i want to include the modification of the equations in order to take into account a variable density and viscosity. Here, no surface tension will be included. The diffuse interface will be at phi = 0.5.
I want to test the case of Reyleigh-Taylor instability

The Attempt at a Solution

Sorry everyone, I am not in this domain, I am doing a thesis in combustion but i need this program in my works.
Thank you in advance for your help, i will appreciate that.

 

[Chestermiller]

Did you really expect anyone to be able to figure out what your coding is all about without any real detailed explanation?

 

[jim mcnamara]

Some points to help you get help:
Please use code tags, and comment your code. Why? Well, for example python is sensitive to leading spaces. Which is a good thing.
Your question is currently a wall of hard-to-read and syntactically incorrect code.

https://www.quora.com/Why-is-Python-space-sensitive

 

[BvU]

Hello mr 1,

As you can read in the PF guidelines, 'Hello everyone' does not constitute a relevant equation, and 'Sorry everyone' does not count as an attempt at solution. So please do better than that.

And under 1. I expect at least a problem statement and I'm not referring to

modify this program to include an scalar phi corresponding to a liquid fraction

Does that mean there is already a tensor $\phi$ in there and you want it simpler ?

And how about a desscription and a list of variables and their meaning ?

I can put $[$code=python$]$ and $[$\code$]$ around your post, but that - of course - does not restore the missing spaces:

import time
import numpy as np
from Output import *
from PoissonSolver import *
tmps1=time.clock()

# Physical parameters
Re = 10.

# Numerical parameters
CFL = 0.9
dt = 0.001/4
t = 0.0
tmax = 1.5

# Make mesh
from MeshGeneration import *
mesh = Mesh(2.0,1.0,50,25)

# Output parameters
dtAff = 0.01 # s

# Boundary condition
# Cavity
#BC = {"type":{"Left": "Wall", "Right": "Wall","Bottom": "Wall", "Top": "MovingWall"},
# "u":{"Top": 1.0},
# "v":{}}
# Forcing Vortex
# BC = {"type":{"Left": "MovingWall", "Right": "MovingWall","Bottom": "MovingWall", "Top": "MovingWall"},
# "u":{"Top": 1.0, "Bottom": -1.0},
# "v":{"Left": 1.0, "Right": -1.0}}
# Channel
BC = {"type":{"Left": "Inflow", "Right": "Outflow","Bottom": "Wall", "Top": "Wall"},
"u":{"Left": 1.0},
"v":{"Left": 0.0},
"p":{"Right": 0.0}}
#BC = {"type":{"Left": "Wall", "Right": "Wall","Bottom": "Outflow", "Top": "Inflow"},
# "u":{"Top": 0.0},
# "v":{"Top": -1.0}}

# Initial condition
from InitNS import *
u, v, p = InitRest(mesh)

if BC["type"]["Left"] == "Inflow":
u[0,:] = BC["u"]["Left"]
if BC["type"]["Top"] == "Inflow":
v[:,mesh.Ny] = BC["v"]["Top"]

# Define fluxes
from NS_fluxes import GRAD, STRESS, TRANS, DIV

# Create Matrix for Poisson Solvre
A, bp = CreatePoissonMatrix(mesh, BC)

savefigNS(u,v,p,t,mesh)

while t<tmax:
ustar, vstar = np.copy(u), np.copy(v)
utrans, vtrans = TRANS(u,v,mesh,BC)
ustress, vstress = STRESS(u,v,mesh,BC)

ustar[1:-1,:] = u[1:-1,:] + dt*(utrans+1.0/Re*ustress)
vstar[:,1:-1] = v[:,1:-1] + dt*(vtrans+1.0/Re*vstress)
b = DIV(ustar,vstar,mesh)
p = PoissonSolver(A,b,bp,mesh)
#p = PoissonSolver(b,BC,mesh)
uprime, vprime = GRAD(p,mesh)
u, v = ustar-uprime, vstar-vprime

t = t + dt
# Plot
if (abs(t % dtAff) < dt):
savefigNS(u,v,p,t,mesh)
tmps2=time.clock()
print "Temps d'execution = %d\n" %(tmps2-tmps1)