# How to solve a system of six coupled PDE with step function

i'm trying to solve a coupled system of six partial differential equations with NDSolve. The system is as follows but i don't have idea how to use the conditionals. Here $H_{\sigma}$ is a continuous mollifier of the step function satisfying The parameter $\sigma > 0$ fixes the thickness of the transition between $H_\sigma(s) = 1$ and $H_\sigma(s) = 0$

My code is below

H[s_, \[Sigma]_] := Which[s <= 0, 1, s > \[Sigma], 0, True, 1 - s/\[Sigma]]

solucion[t_] =
With[{gaman = 0.746, gamat = 0.97, psi0 = 0.75, psi = 0.655,
psin = 0.6, psit = 0.8, deltat = 0.03, deltan = 0.1, mun = 0.1,
mut = 0.05, nu = 0.000016, pin = 6000000,
pit = 3000000, \[Tau] = 0.005, fo = 0.25, fg = 0.16, betan = 1.3,
betat = 1.2, Kan = 0.1, Kat = 0.3, \[Kappa] = 0.00000734, D = 1.0},

NDSolve[{
fin'[t] == oxi[t]*gaman*fin[t]*H[s, \[Sigma]] (psi - psin) - ((1 -
oxi[t]))*deltan*fin[t],

fit'[t] == glu[t]*gamat*fit[t]*H[s, \[Sigma]] (psi - psit) - ((1 -
glu[t]))*deltat*fit[t],

mat'[t] == mun*fin[t] + mun*fit[t] - nu*enz[t]*mat[t],

enz'[t] == pin*fin[t] + pit*fit[t],

oxi'[t] == -betan*fin[t]*oxi[t] + fo,

glu'[t] == -betat*fit[t]*glu[t] + fg,

fin == 0.45, fit == 0.005, mat == 0.2, enz == 0.2,
oxi == 0.25, glu == 0.16}, {fin[t], fit[t], mat[t], enz[t],
oxi[t], glu[t]}, {t, 0, 80}] /. {s -> 0, \[Sigma] -> 1}][[1, All, 2]]

Plot[solucion[t], {t, 0, 80}, PlotRange -> {{0, 80}, {0, .55}}]


but the above code does not correctly plot the expected results.

My explicit question is how to properly use conditionals inside NDSolve for $H_{\sigma}(s)$ for those three cases??

Is it possible? Could you please help me?

Any suggestions are very appreciated

*The image below corresponds to the same eqs plotted in Python. • Are you sure the definition of $H_{\sigma}(s)$ is correct? Its value is $1$ as long as $s\neq0$? Just try Plot3D[Evaluate@H[x, y], {x, -5, 5}, {y, -5, 5}] and you'll know what I mean. – xzczd Oct 25 '17 at 3:04
• check the definition of the step function and I was wrong. I uploaded the new image of the definition of H thx for the observation – Rodrigo López Oct 25 '17 at 4:37
• What's the expected result? I think NDSolve has solved the given equation system correctly. – xzczd Oct 25 '17 at 5:21
• I have a code made in Python and the results plotted are different. e.g. the first eq in the system "fin'[t]" plotted in python has the expected behavior as shown in the image (i uploaded a new image plotted in python). I think the error is in properly using the step function. In addition there is a message "NDSolve :: ndnum: Encountered non-numerical value for a derivative at t == 0.", however, I have no idea why it shows it. – Rodrigo López Oct 25 '17 at 19:24
• 1. The warning is simply because you don't place ReplaceAll (/.) in the correct place, but it doesn't cause mistake in the final result. If you don't want to see the warning, place the /. {s -> 0, \[Sigma] -> 1} inside NDSolve. 2. The symbol in your code is quite different from those in the picture, coding equations in this way is not a good idea because this makes error checking painful and inefficient. 3. The 4th equation is apparently different from the one in the picture. 4. You'd better add the Python code in your question if it's not too long. – xzczd Oct 26 '17 at 3:57