usnistgov / fipy

FiPy is a Finite Volume PDE solver written in Python
http://pages.nist.gov/fipy/en/latest
Other
504 stars 148 forks source link

The coefficient must be a vector value. #940

Closed ghost closed 1 year ago

ghost commented 1 year ago

I've been trying to simulate a 2D, two phase flow system by using Navier-Stokes equations but, I keep tackling with this error 'The coefficient must be a vector value." I've found inspiration in https://github.com/usnistgov/fipy/blob/master/examples/reactiveWetting/liquidVapor2D.py. I've managed to solve this error but I couldn't! Does anybody know how I can debug the error? The error is: VectorCoeffError Traceback (most recent call last) ~\AppData\Local\Temp\ipykernel_12276\1560868714.py in 38 #Define the Navier-Stokes equations 39 ---> 40 equ1 = TransientTerm(var=velocityX) == - ConvectionTerm(coeff=density.faceValue, var=velocityX) + pressure.grad + (waterDensity g) 41 equ2 = TransientTerm(var=velocityY) == - ConvectionTerm(coeff=density.faceValue, var=velocityY) + pressure.grad + (liquidDensity g) 42 equ3 = ConvectionTerm(var=velocityX)

C:\ProgramData\Anaconda3\lib\site-packages\fipy\terms\abstractConvectionTerm.py in init(self, coeff, var) 75 76 if isinstance(coeff, _MeshVariable) and coeff.rank < 1: ---> 77 raise VectorCoeffError 78 79 if isinstance(coeff, CellVariable):

VectorCoeffError: The coefficient must be a vector value.

The code is: `#Import fipy from fipy import CellVariable, Grid2D, TransientTerm, VanLeerConvectionTerm, DiffusionTerm, ImplicitSourceTerm, ConvectionTerm, CentralDifferenceConvectionTerm, Viewer

Set up the domain

L1 = 20.0 #cm L2 = 5.0 #cm

x = 10.0 y = 2.0

dx = L1 / x dy = L2 / y

mesh = Grid2D(nx=x, ny=y, dx=dx, dy=dy) x, y = mesh.cellCenters

Define the initial & physical parameters

g = 981.0 #cm/s2 waterDensity = 1.0 #g/cm3 liquidDensity = 0.85 #g/cm3

Define the fluid velocity

velocityVector = CellVariable(mesh=mesh, name='Velocity', rank=1) velocityX = CellVariable(mesh=mesh, hasOld=True, name='Velocitx') velocityY = CellVariable(mesh=mesh, hasOld=True, name='Velocity')

Define the density & pressure

density = CellVariable(mesh=mesh, hasOld=True, name='Density') pressure = CellVariable(mesh=mesh, hasOld=True, name='Pressure')

Define boundary condition

velocityX.constrain(1, mesh.facesLeft & mesh.facesRight) velocityY.constrain(0, mesh.facesTop & mesh.facesBottom)

Define the Navier-Stokes equations

equ1 = TransientTerm(var=velocityX) == - ConvectionTerm(coeff=density.faceValue, var=velocityX) + pressure.grad + (waterDensity g) equ2 = TransientTerm(var=velocityY) == - ConvectionTerm(coeff=density.faceValue, var=velocityY) + pressure.grad + (liquidDensity g) equ3 = ConvectionTerm(var=velocityX) equ4 = ConvectionTerm(var=velocityY)

Solve the equations

eqns = equ1 & equ2 & equ3 & equ4 timeStepDuration = 0.1 steps = 10 for step in range(steps): eqns.solve(dt=timeStepDuration)

Visualize the results

viewer = Viewer(vars=(density, velocityX), datamin=0.0, datamax=1.0) viewer.plot()`

Acid_Model_13.zip