Skip to main content
5-Regular Member
February 2, 2026
Question

ODE's error with if else function

  • February 2, 2026
  • 1 reply
  • 226 views

I am solving a system of 1D steady ODEs in Mathcad using RK4. The system works correctly until I add a heat transfer term representing two-way coupling between a gas phase and injected droplets. When this term is included, the solver no longer outputs results and reports that the gradient matrix must be scalar or a matrix.

The gas stagnation temperature gradient depends on Reynolds number, Nusselt number and a convective heat transfer coefficient that is itself a function of the ODE variables, which vary along z (length of nozzle)

All other gradients (Mach, pressure, velocity, etc.) are functions of the same state vector and work correctly when the heat term is removed. I suspect I am incorrectly defining a function within the gradient matrix or using If else wrong.

Any help would be appreciated, as I'm going round in circles.

below is screenshot of part of code:

MC_14395780_3-1770047934724.png

 

MC_14395780_0-1770047771623.pngMC_14395780_1-1770047833792.pngMC_14395780_2-1770047863631.png

 

Thanks very much

 

1 reply

25-Diamond I
February 2, 2026

How about posting the worksheet itself instead of just pictures? And clarifying exactly which added term made the system fail.

 

What I noticed is that in your definitions of some grad... functions you divide by x7(y1) and in the vector x you provide x7 is zero.

Furthermore it looks like x8(T01) never is used, but I may have overlooked it in the pictures.

5-Regular Member
February 2, 2026

Hi, 

I have now added the file above. I have removed x8 (T01) from the vector, and changed y1 (x7) to 0.1.

I think it may be something to do with my if else statement but I can't spot the mistake.

25-Diamond I
February 2, 2026

You can use Prime's error trace facility to trace back to the position the error shows up first.

 

In your case it traces back to the call of function F.drag.

You defined F.drag as a function with two arguments

Werner_E_0-1770059833561.png

but you call it with just one argument

Werner_E_1-1770059873172.png

The error message tells you exactly that 😉

 

Same problem with function C.d

Werner_E_2-1770060019135.png

Unfortunately fixing this

Werner_E_3-1770060342680.png

don't help for rkfixed.

It looks like the results do not depend much on z. The results are the same for z=0 and z=L.e=0.5.
Are the values expected to that high? (10^11)

 

Which was omitted when the system finally worked OK? What was added causing the problem?