Hello,
Being a (relatively) new DuMux user, I again ran into a problem that I
can not understand.
Long story short, I once had a working TwoPNC simulation (two
incompressible liquids, one-component and two-component, with diffusion
set to zero). Then I tried to add temperature to it, using TwoPNCNI, and
failed. Then I simplified the setup as much as I can and it still failed.
What I have now is: 1x1 grid, Neumann conditions on all faces. Zero
liquid fluxes and a nonzero heat flux are prescribed though the left
boundary. On the right boundary, a zero heat flux is prescribed. Outflux
of liquid can exist at the right boundary, because of the fixed pressure
there, but in fact the fixed boundary pressure is equal to the (initial)
pressure inside.
So what I naturally expect in this example is that there should be no
liquid flux, while the temperature should rise. What I get instead is
some numerical problem. Here is the output with some debug printing
which I added into NewtonSolver::solveLinearSystem_:
Solve: M deltax^k = r
jacobian [n=1,m=1,rowdim=4,coldim=4]
row 0 0.000000000000000e+00 -1.700000211712904e+12
-3.400002242415212e+11 0.000000000000000e+00
row 1 1.166881779343447e-07 7.431428154328993e+10
0.000000000000000e+00 0.000000000000000e+00
row 2 0.000000000000000e+00 0.000000000000000e+00
2.692656567250088e+13 0.000000000000000e+00
row 3 -3.075568914721009e+07 1.403093338012695e+15
0.000000000000000e+00 2.398149408799923e+14
residual = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 8.000000000000000e+02]
deltaU = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 3.335905582297870e-12]
result = 1
Newton iteration 1 done, residual = 6.066e+00
Solve: M deltax^k = r
jacobian [n=1,m=1,rowdim=4,coldim=4]
row 0 0.000000000000000e+00 -1.700000211712904e+12
-3.400002242415212e+11 0.000000000000000e+00
row 1 1.166881779343447e-07 7.431428154328993e+10
0.000000000000000e+00 0.000000000000000e+00
row 2 0.000000000000000e+00 0.000000000000000e+00
2.692656567250088e+13 0.000000000000000e+00
row 3 -3.075568914721009e+07 1.403093338012695e+15
0.000000000000000e+00 2.398146891623329e+14
residual = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 6.066131591796875e+00]
deltaU = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 2.529507935058412e-14]
result = 1
Newton iteration 2 done, residual = 6.066e+00
Solve: M deltax^k = r
jacobian [n=1,m=1,rowdim=4,coldim=4]
row 0 0.000000000000000e+00 -1.700000211712904e+12
-3.400002242415212e+11 0.000000000000000e+00
row 1 1.166881779343447e-07 7.431428154328993e+10
0.000000000000000e+00 0.000000000000000e+00
row 2 0.000000000000000e+00 0.000000000000000e+00
2.692656567250088e+13 0.000000000000000e+00
row 3 -3.075568914721009e+07 1.403093338012695e+15
0.000000000000000e+00 2.398146891623329e+14
residual = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 6.066131591796875e+00]
deltaU = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 2.529507935058412e-14]
Newton iteration 3 done, residual = 6.066e+00
[...]
Newton iteration 20 done, residual = 6.066e+00
Newton solver did not converge with dt = 0.001 seconds. Retrying with
time step of 0.0005 seconds
Solve: M deltax^k = r
jacobian [n=1,m=1,rowdim=4,coldim=4]
row 0 0.000000000000000e+00 -3.400000423425809e+12
-6.800004484830424e+11 0.000000000000000e+00
row 1 1.166881779343447e-07 1.486285630865799e+11
0.000000000000000e+00 0.000000000000000e+00
row 2 0.000000000000000e+00 0.000000000000000e+00
5.385313134500177e+13 0.000000000000000e+00
row 3 -6.151137829442018e+07 2.806186676025390e+15
0.000000000000000e+00 4.796298817599846e+14
residual = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 8.000000000000000e+02]
Newton: Caught exception: "FMatrixError
[luDecomposition:/home/dpavlov/DUMUX2/DUMUX/dune-common/dune/common/densematrix.hh:909]:
matrix is singular"
[...]
Newton solver did not converge with dt = 1.95313e-06 seconds. Retrying
with time step of 9.76563e-07 seconds
Solve: M deltax^k = r
jacobian [n=1,m=1,rowdim=4,coldim=4]
row 0 0.000000000000000e+00 -1.740800216794014e+15
-3.481602296233177e+14 0.000000000000000e+00
row 1 1.166881779343447e-07 7.609782430032889e+13
0.000000000000000e+00 0.000000000000000e+00
row 2 0.000000000000000e+00 0.000000000000000e+00
2.757280324864090e+16 0.000000000000000e+00
row 3 -3.149382568674313e+10 1.436767578125000e+18
0.000000000000000e+00 2.455704994611121e+17
residual = [0.000000000000000e+00 0.000000000000000e+00
0.000000000000000e+00 8.000000000000000e+02]
Newton: Caught exception: "FMatrixError
[luDecomposition:/home/dpavlov/DUMUX2/DUMUX/dune-common/dune/common/densematrix.hh:909]:
matrix is singular"
Dune reported error: NumericalProblem
[solve:/home/dpavlov/DUMUX2/DUMUX/dumux/dumux/nonlinear/newtonsolver.hh:254]:
Newton solver didn't converge after 10 time-step divisions. dt=9.76563e-07
I see that the matrix is not in great condition, so to speak, but what
can I do about it? Is my setup wrong? Does the DiffMethod::numeric
assembler work bad? Should I define more stuff about thermal properties
of my system?
Right now, I have liquidHeatCapaticity(), liquidThermalConductivity(),
liquidEnthalpy() defined for my liquids, and also I set
ThermalConductivityModel property to ThermalConductivityJohansen.
Finally, I set Component.SolidDensity,
Component.SolidThermalConductivity, and Component.SolidHeatCapacity
parameters.
Best regards,
Dmitry
_______________________________________________
Dumux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
