> > Hi, > > I tried to solve a problem with a Robin Boundary condition in fipy. The > physical problem is very simple, but I could not get the numerical solution > to match the analytical solution. > > [image: Inline image 2] > > Is there something wrong with the way I set up the boundary condition? > I've copied my script below and bolded the boundary condition lines. Thanks > for your help. > > Best regards, > Rose > ---------------------------------------------- > > from fipy import * > > L = 1. > > nx = 100 > > dx = L/nx > > mesh = Grid1D(nx=nx,dx=dx) > > zc = mesh.cellCenters[0] > > Dm = 8e-6 # moleculer diffusion in water column > > Dmp = 5e-6 # moleculer diffusion in porewater (corrected for tortuosity) > > Dbw = 1e-5 # bioturbation in porewater > > D = Dmp+Dbw > > > Cw=1e-3 > > Cl=1. > > > time = Variable() > > beta = 1e-3 > > #*abs(numerix.sin(w*time)) > > > # steadystate solution > > A = beta*(Cl-Cw)/(D+beta*L) > > B = Cl-A*L > > Canal = CellVariable(name = "analytical solution steadystate", mesh=mesh, > value = A*zc+B) > > > C = CellVariable(name="concentration", mesh = mesh,value=0.) > > *C.faceGrad.constrain(beta*(C.faceValue-Cw)/D, mesh.facesLeft)* > > *C.constrain(Cl,mesh.facesRight)* > > eqI =TransientTerm()==ImplicitDiffusionTerm(coeff=D) > > > viewer = Viewer(vars= (Canal,C)) > > dt = 864 > > while time()<86400: > > time.setValue(time() + dt) > > eqI.solve(var = C,dt=dt) > > viewer.plot() >
_______________________________________________ fipy mailing list [email protected] http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
