0

With FiPy (Python library) I want to solve the coupled of pdes system shown below. The code below works, but doesn't give the correct solution.

p, q, r, s = 2, 1, 2, 0
Du, Dv = 0.0004, 0.04
mesh = Grid1D(nx=500, Lx=5.)

U = CellVariable(name = 'U',mesh=mesh, hasOld=True, value=1.)
V = CellVariable(name = 'V',mesh=mesh, hasOld=True, value=2.)

eqU = TransientTerm(var=U) == DiffusionTerm(coeff=Du, var=U) - U + (U**p)/(V**q)
eqV = TransientTerm(var=V) == DiffusionTerm(coeff=Dv, var=V) - V + (U**r)/(V**s)

viewerV = Viewer(vars=V)
viewerU = Viewer(vars=U)

timeStepDuration = 0.1
steps = 1000

eqn = eqU & eqV

for t in range(500): 
    U.updateOld()
    V.updateOld()
    eqn.solve(dt=1.e-3)
viewerV.plot()
viewerU.plot()

In line with some code I saw I also tried to replace the CellVariables with ImplicitSourceTerm, but this gives an error, since taking a ImplicitSourceTerm to the power p can't be done in this way. However I can't find any documentation on how it should be done:

eqU = TransientTerm(var=U) == DiffusionTerm(coeff=Du, var=U) - ImplicitSourceTerm(var=U) + ((ImplicitSourceTerm(var=U)**p)/(ImplicitSourceTerm(var=V)**q))
eqV = TransientTerm(var=V) == DiffusionTerm(coeff=Dv, var=V) - ImplicitSourceTerm(var=V) + ((ImplicitSourceTerm(var=U)**r)/(ImplicitSourceTerm(var=V)**s))

Trying to raise U.value**p also gives errors.

1 Answers1

0

There's nothing intrinsically wrong with your first solution, other than that the explicit sources may not converge as rapidly as implicit sources would.

An ImplicitSourceTerm represents some coefficient times the variable that's being solved for (the default coefficient is 0, which is useless, so you must explicitly state coeff=1.).

Your terms could be expressed as

eqU = TransientTerm(var=U) == DiffusionTerm(coeff=Du, var=U) - ImplicitSourceTerm(var=U, coeff=1.) + ImplicitSourceTerm(var=U, coeff=(U**(p-1))/V**q)
eqV = TransientTerm(var=V) == DiffusionTerm(coeff=Dv, var=V) - ImplicitSourceTerm(var=V, coeff=1.) + ImplicitSourceTerm(var=U, coeff=(U**(r-1))/V**s)

You should then use sweeps.

jeguyer
  • 2,379
  • 1
  • 11
  • 15
  • Thanks for the responds! I still have trouble getting the right results. When I set the coeff to U**(p-1), I get the error: TermMultiplyError: Must multiply terms by int or float. Using: f = (U**(r-1))/(V**s)-1. g = (U**(p-1))/(V**q) Runs, but doesn't give correct results. I have been playing around with various ways to implement this term but all failed. The system either gives U=V=1 as a solution for the whole grid or it goes to infinity for the whole grid. Where it should show a 'spatial' sinusoid over the grid. – Viktor van der Valk Mar 15 '18 at 21:59
  • Also changing U to ImplicitSourceTerm(var=U, coeff=1.) gives completely different solutions. – Viktor van der Valk Mar 15 '18 at 22:04
  • Sorry about that. The denominator `V**q` or `V**s` should be part of `coeff`, rather than dividing the `ImplicitSourceTerm`, which isn't allowed. – jeguyer Mar 19 '18 at 15:56
  • I am not surprised that changing `U` to `ImplicitSourceTerm(var=U, coeff=1.)` changes the results. Are you sweeping? – jeguyer Mar 19 '18 at 15:57