• Problem setup
• Workaround

C code generation with intermediate complex expressions

Christopher D. Morgan, P.E., Ph.D.

Sometimes, a calculation requires intermediate expressions containing complex numbers, even though the final result is a real number. This can create some problems when generating Fortran or C routines. In this application, I demonstratate a workaround to this difficulty.

> restart;

Problem setup

Begin by defining a real number, computed by multiplying a complex number with its conjugate.
Notice that xx is real.
Another way to build x_cc --required if a and/or b are complex-- would be to use the "conjugate" function. That level of generality is not required for code generation since a and b must be real.
x := cos( a+I*b)^5:
x_cc:= conjugate( x );
xx:= x*x_cc;

Convert xx to a procedure.

> x_function:=codegen[makeproc]( xx, [a, b] );

Convert x_function procedure to C. This doesn't work because the expression contains the complex number i, even though the final result is purely real.
codegen[C] ( x_function );

Error, (in codegen/C) unable to translate I

Workaround

Use evalc to separate the expression into its real and imaginary parts. The imaginary parts sum to zero, but Maple has not yet performed this simplification.

> x1:= evalc( xx);

Simplify the above expression so that the imaginary terms cancel each other. This is often --but not always-- needed.

> x1:=simplify( x1 );

After using evalc and simplify, the rest of the code generation is just as expected.

> x_function:=codegen[makeproc]( x1, [a, b] );

> codegen[C] ( x_function );

#include <math.h>

double x_function(a,b)

double a;

double b;

{

  {

    return(-1.0-20.0*pow(cos(a),2.0)*pow(cosh(b),2.0)-30.0*pow(cos(a),4.0)*pow(

cosh(b),4.0)+5.0*pow(cos(a),2.0)+5.0*pow(cosh(b),2.0)-10.0*pow(cos(a),4.0)-10.0

*pow(cosh(b),4.0)+pow(cosh(b),10.0)+pow(cos(a),10.0)-5.0*pow(cos(a),8.0)-5.0*

pow(cosh(b),8.0)+10.0*pow(cos(a),6.0)+10.0*pow(cosh(b),6.0)+30.0*pow(cos(a),2.0

)*pow(cosh(b),4.0)+30.0*pow(cos(a),4.0)*pow(cosh(b),2.0)-20.0*pow(cos(a),2.0)*

pow(cosh(b),6.0)+10.0*pow(cos(a),4.0)*pow(cosh(b),6.0)+10.0*pow(cos(a),6.0)*pow

(cosh(b),4.0)-20.0*pow(cos(a),6.0)*pow(cosh(b),2.0)+5.0*pow(cos(a),8.0)*pow(

cosh(b),2.0)+5.0*pow(cos(a),2.0)*pow(cosh(b),8.0));

  }

}