Nonlinear Equations: How to Solve Deviations from Textbook Solutions

[DrScanlon]

First off, thanks to the people who started this forum, and I apologize for this post being so long-winded.

I'm a retired engineer dabbling in engine simulation. The text I'm using as the foundation for the gas dynamics theory contains several nonlinear equation sets for each simulated scenario, followed by the answers for the unknown variables.

Using one particular set of 4 equations, I use the Newton-Raphson (N-R) method, as recommended by the author, and obtain the solution for the 4 unknowns. The problem is that my solution (the variable values) deviates substantially from the solution in the text even if I use the author's solved values--or something very near them in value--as my initial first guess. I also see this problem when using a commercial nonlinear solver. I thought a good test would be to plug the text solution values into each of the equations, expecting to see a "0" result for each, but instead, some results were far from zero. All of the coefficients match those in the text and I have verified that the equations are correct.

I trust the answers in the text; however, I suspect that something other than the plain vanilla N-R method is used, or that there is some "best fit" solution that doesn't require the equations to resolve to "0". I have sufficient intuitive knowledge of the system model that I could use some form of bracketing on the first guess for the unknowns, but I have no clue in this area.

Can anyone recommend a next step in trying to resolve this?

Frank

 

[TheoMcCloskey]

DrScanlon - don't be too quick to trust the book's answers - they are not infallible - especially when you find (significant) non-zero values when you pluged those "solutions" into the original equations.

I unfortunately don't have any free time at the momment but I would be interested in reviewing the problem and so-called solutions. Can you post the equations and the author's proposed solutions?

 

[DrScanlon]

The text numerical results could be in slight error by a few percent, but they make physical sense. Plugging the answers back into the N-R yields values that are not physically possible, considering the input.

I'll post the equations, as well as the input values and supporting explicit equations. I'm assuming that since you may want to cut/paste for testing, that it might be best to avoid Tex format for the equations.

Thanks for the reply.

 

[TheoMcCloskey]

Thanks - Yes, I am curious and would like to review, if you don't mind, to see if I see anything obvious as to what you're experiencing.

Sometimes these equations can be very unstable near the roots and may tend to drift off to spurious "solutions" away from the "solution of interest" but I'm somewhat doubtful that this is the case here. But, let's take a look anyway.

Either format is fine, please just make sure it accurately represents the problem.

Thanks again

 

[DrScanlon]

Theo-- The expressions that I input to the N-R solver, as well as the supporting data/equations are below. Some vars may be defined, but not used. I pulled this from my code, so if there are errors, please let me know and I'll fix them. Again, thanks for analyzing this.

Frank

// ! EQUATIONS/DATA -- BEGIN !

Expressions ( all equal "0")

eq0:
((p0*X1^2)/(R*T1))*(Xt^G5)*At*ct
-((p0*gamma)/(a02^2))*(Xi2+Xr2-1)^G5*A2*G5*a02*(Xi2-Xr2)

eq1:
G5*(a01*X1)^2-((G5*a02*(Xi2-Xr2))^2 + G5*a02^2*(Xi2+Xr2-1)^2)

eq2:
G5*((a01*X1)^2 - (a01*Xt)^2) - ct^2

eq3:
p0*At*(Xt -(Xi2+Xr2-1)^G7) +(((p0*gamma)/(a02^2))
* A2*((Xi2+Xr2-1)^G5)*G5*a02*(Xi2-Xr2))*(ct-G5*a02*(Xi2-Xr2))

Unknown Variable Initialization (first guesses)
(NOTE: These are also the solution values from the text)

Xr2 = 1.0577
Xt = 1.059
a02 = 522.1
ct = 262.9


Input Data

P1 = 1.8
T1 = 500
dt = 25.0
d2 = 30.0
cd = .75
Pi2 = 1.1


Explicit Equations and Constants
(NOTE: Assumes Input Data is defined)
R = 300
gamma = 1.36
p0 = 101325
g5 = 2.0/(gamma - 1)
g7 = (2*gamma)/(gamma - 1)
g15 = 1.0/g5
g17 = 1.0/g7
p0 = 101325.0
Xi2 = Pi2^g17
X1 = P1^g17
T01 = (T1+273.0) /X1^2
a01 = sqrt(gamma * R * T01)
rho_01 = p0 / (R * T01)
At = cd*dt^2* 0.7854
A2 =dt^2*0.7854

// ! EQUATIONS/DATA -- END !

 

[DrScanlon]

Thanks - Yes, I am curious and would like to review, if you don't mind, to see if I see anything obvious as to what you're experiencing.

Sometimes these equations can be very unstable near the roots and may tend to drift off to spurious "solutions" away from the "solution of interest" but I'm somewhat doubtful that this is the case here. But, let's take a look anyway.

Either format is fine, please just make sure it accurately represents the problem.

Thanks again

Theo-- I've added a correction to the equations as follows:

Correction #1:

eq3 - expression is changed to read as follows:

p0*(Xt^G7 -(Xi2+Xr2-1)^G7)- ((p0*gamma)/(a02^2))*(Xi2+Xr2-1)^G5*G5*a02*(Xi2-Xr2)
*(ct-(G5*a02*(Xi2-Xr2)))

This is straight from the text. It is the momentum equation which usually includes an Area (A) term. I did a quick derivation and could not get the "A" terms to drop out; hence, I was using the one with the A terms in my code. With that said, I believe it's best to start with the originals from the text as a clean starting point. I'm going to do a more thorough review of this equation in an attempt so see where the area term(s) drop out in the derivation.

Frank