Summary

Newton's method is a numerical method for finding the zeros of non-linear equations. There are many applications that arise in CFD that can benefit from this technique. This article provides two such examples and outlines how to use simulation operations in STAR-CCM+ to automate the processes. Two techniques are discussed, with the use of the adjoint solver and without.

Details

Attachments:

CRM_Airfoil_Adj_NM_TargetLift.sim (41 MB) UCAV_FindMdotBC_NM_Coarse.7z (1 GB)

Using Newton's Method with Simulation Operations in STAR-CCM+

Newton's method is a well-known technique for finding the zeros of non-linear equations.  Many situations arise in CFD calculations that can be framed as an exercise in finding the zeros of such equations and Newton’s method can be used to efficiently solve these problems.

Take for example the process of identifying the correct intake and exhaust mass flow boundary conditions for external aerodynamics simulations of aircraft with air breathing engines.  At steady state flight conditions the drag and thrust forces should balance.
 

 

At fixed Mach and altitude, the engine exhaust mass flow rate  is used to balance these forces:
 

 

We therefore seek a value of  such that  which is the same as looking for the zeros of the function  for which we can use Newton’s method.

Newton's method for a function of a single variable can be derived from a Taylor series expansion:
 

.
 

Truncating the higher order terms and setting  we can then solve for
 

 

This equation is Newton's method, and it is solved iteratively until the value of  is within some user specified tolerance of 0. The subject of this article is to implement this method in a STAR-CCM+ simulation.  In some situations, the value of  can be found using the adjoint solver but when this is not the case we require an approximation.  Using backward differencing we obtain the following equation which can be implemented in STAR-CCM+ using simulation operations.
 

 

When the backwards difference approximation is used in place of the derivative we are now using the Secant Method instead of Newton's method.  Two initial guesses  and  are required. Our implementation will rely on the use of simulation parameters, expression reports, and simulation operations.  All of these objects are tagged as "Newtons Method" in the attached example sim files and an associated custom tree is available in the tools section to make it easy to locate all the relevant objects.

We will deploy the method to determine the appropriate engine exhaust mass flow boundary condition of a generic drone external aerodynamics simulation.
 

 

In the animations below we can see how the Secant method progresses from our initial two guesses of  and . 
 

 

With each subsequent iteration the value of  gets closer to 0. After 4 iterations it is less than 1% of the drag value. In all 6 values for  were tried (2 initial guesses and 4 from the Secant method) requiring less than 800 iterations to complete the process.  A description of the steps in the simulation operation is given below.  This is taken from the attached example simulation UCAV_FindMdotBC_NM_Coarse.7z

 

 

For the first iteration after a new value of  has been  selected the values of  and   are the same which can result in a divide by zero.  While this won’t have an impact on the results of the simulation it can cause the simulation operation manager to detect an error and terminate when it shouldn’t.  This can be avoided by using an expression to check that these two values are not equal before using them to compute the new value of .  For example in the report "dmdot_dF_net" used to calculate :
 

abs(${F_net} - ${F_net_0}) < 0.0000001 ? 1 : (${mdot_out_1} - ${mdot_out_0}) / (${F_net} - ${F_net_0})

 

Newton's method when using the adjoint solver

If you can use the adjoint solver in STAR-CCM+ to compute the derivative  then the implementation can be simplified.  As an example, suppose you have an airfoil, and you would like to know at what angle of attack  the lift coefficient .  Here our function  is the difference between the current  and the desired value of 1.025 or

 

 

Newton's method in this case becomes
 

 

Where the derivative  can be computed using the adjoint solver (please see article Example of using adjoint to calculate the gradient of lift with respect to angle-of-attack, or attached sim file CRM_Airfoil_Adj_NM_TargetLift.sim ), so the backwards difference approximation is not necessary and only a single initial guess is required. In this case the simulation operations can be simplified to the following:
 

After the initial guess of  two additional iterations were required to identify that  when .
 

 

Newton's method is expected to converge as long as the initial guess(es) are close enough to the zero of interest.  In most CFD applications this is not a challenging requirement, but in cases where the function in question has many zeros or has many local minima or maxima a good initial guess(es) may be necessary.