This will be the last of the series on solving non-linear equations (for now). Up until now all the examples have had two unknown values, and two target values. This can be extended by making three changes to the code:
- Set up an nxn matrix of the function slopes with respect to each of the unknown values (the function Jacobian)
- Calculate the first estimate of the target values, using an estimated value for each of the unknowns.
- Form and solve a series of simultaneous equations to find a better estimate of the unknowns.
In Excel VBA the solution to the equations can be found using the Worksheetfunction.MInverse and MMult functions, as shown in the code below:
LoopNum = LoopNum + 1 ' Evaluate function at guessed values and small increments. For i = 1 To NumVar Var2A(i, 1) = Var1A(i, 1) * VarFact VarDiffA(i, 1) = Var1A(i, 1) * (VarFact - 1) Next i Res1 = Application.Run(Func, Var1A, X) For j = 1 To NumVar ResDiff(j, 1) = Target(j, 1) - Res1(j, 1) Next j Temp = Var1A(1, 1) Var1A(1, 1) = Var2A(1, 1) For i = 1 To NumVar Res2 = Application.Run(Func, Var1A, X) For j = 1 To NumVar SlopeA(j, i) = (Res2(j, 1) - Res1(j, 1)) / VarDiffA(i, 1) Next j Var1A(i, 1) = Temp If i < NumVar Then Temp = Var1A(i + 1, 1) Var1A(i + 1, 1) = Var2A(i + 1, 1) End If Next i ' Solve SlopeA InvA = WorksheetFunction.MInverse(SlopeA) ResA = WorksheetFunction.MMult(InvA, ResDiff) ErrSum = 0 For i = 1 To NumVar ErrA(i, 1) = (Target(i, 1) - Res1(i, 1)) ErrSum = ErrSum + Abs(ErrA(i, 1)) Var1A(i, 1) = Var1A(i, 1) + ResA(i, 1) Next i
The User Defined Function (UDF) MSolve has been modified as shown above, and the earlier version renamed MSolve2. A spreadsheet including examples of both functions, and full open source code, can be downloaded from:
Input and output from both functions is shown in the screenshot below:
MSolve2 works as before, in the example returning the depth of neutral axis and compression face strain in a reinforced concrete section, for any given axial load and bending moment.
The MSolve example adds one more unknown (the diameter of the tension reinforcement), and another target, the tension reinforcement stress.
The Python Scipy package also contains a number of routines for solving problems of this type. In order to access these from Excel we need:
- A python routine for the function to be solved.
- An interface allowing this routine and the Scipy functions to be accessed from Excel.
I have updated the Eval-PyIntegration spreadsheet to perform this task, and this can be downloaded from:
The spreadsheet uses ExcelPython for communication with Python, and the download file includes all necessary files (other than Excel, Python and Scipy).
The example below shows the xl_SolveF function used to find the depth of neutral axis, compression face strain, and tension bar diameter required for a specified axial force, bending moment and tension steel stress. In this example the Python function being solved (RCForceMS) is the simplified version, ignoring concrete in tension.
The download also includes Python versions of the functions CurveatMA and CurveatMAS, with input as shown below. The functions make use of the SciPy optimize.root function, which has many optional arguments, detailed in the SciPy manual. In the spreadsheet the only options included are for the solver type, and tolerance for termination.
The screenshot below shows the results of CurveatMA for an axial load of zero and bending moment of zero. The returned values of depth of neutral axis (DNA) and compression face strain (Epsc) have been checked using the function SectForceMV, which finds that the axial force and bending moment are as expected.
The xl_SolveF function can also be called directly, using a function that returns the difference of the function values from the target values, in this case SectForceMDiff.
Finally, CurveatMAS adjusts the neutral axis depth, compression face strain, and tension steel area such that the axial force, bending moment and section curvature match the target values. Once again, xl_SolveF may be used with SectForceMdiff to return the same 3 values.
Note that the steel force and tensile strain found by these functions are average values, including tension stiffening effects. The results are intended for use in deflection calculations, and should not be used in strength or stress limit calculations.