# include # include # include # include # include # include "zero_muller.h" int main ( ); void test01 ( ); void test02 ( ); void test03 ( ); double complex func01 ( double complex x ); double complex func02 ( double complex x ); double complex func03 ( double complex x ); void timestamp ( ); /******************************************************************************/ int main ( ) /******************************************************************************/ /* Purpose: zero_muller_test() tests zero_muller(). Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt */ { timestamp ( ); printf ( "\n" ); printf ( "zero_muller_test():\n" ); printf ( " C version\n" ); printf ( " Test zero_muller(), which uses Muller's method,\n" ); printf ( " with complex arithmetic, to solve a nonlinear equation.\n" ); test01 ( ); test02 ( ); test03 ( ); /* Terminate. */ printf ( "\n" ); printf ( "zero_muller_test():\n" ); printf ( " Normal end of execution.\n" ); printf ( "\n" ); timestamp ( ); return ( 0 ); } /******************************************************************************/ void test01 ( ) /******************************************************************************/ /* Purpose: test01() tests zero_muller() on F(X) = X*X+9. Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt */ { double fatol; double complex fxnew; int itmax; double complex x1; double complex x2; double complex x3; double complex xnew; double xatol; double xrtol; printf ( "\n" ); printf ( "test01():\n" ); printf ( " Demonstrate zero_muller() on F(X) = X*X+9.\n" ); fatol = 1.0E-05; itmax = 10; x1 = 1.0 + 0.0 * I; x2 = 0.0 + 1.0 * I; x3 = 0.5 + 0.5 * I; xatol = 1.0E-05; xrtol = 1.0E-05; zero_muller ( func01, fatol, itmax, x1, x2, x3, xatol, xrtol, &xnew, &fxnew ); printf ( "\n" ); printf ( " X = %20.10f, %20.10f\n", creal ( xnew ), cimag ( xnew ) ); printf ( "\n" ); printf ( " with function value F(X):\n" ); printf ( "\n" ); printf ( " FX = %20.10f, %20.10f\n", creal ( fxnew ), cimag ( fxnew ) ); printf ( " ||FX|| = %20.10f\n", cabs ( fxnew ) ); return; } /******************************************************************************/ void test02 ( ) /******************************************************************************/ /* Purpose: test02() tests zero_muller() on F(X) = (X*X+4) * (X-10) * (X+20) Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt */ { double fatol; double complex fxnew; int itmax; double complex x1; double complex x2; double complex x3; double complex xnew; double xatol; double xrtol; printf ( "\n" ); printf ( "test02():\n" ); printf ( " Demonstrate zero_muller() on F(X) = (X*X+4)*(X-10)*(X+20).\n" ); fatol = 1.0E-05; itmax = 10; x1 = 1.0 + 0.0 * I; x2 = 0.0 + 1.0 * I; x3 = 0.5 + 0.5 * I; xatol = 1.0E-05; xrtol = 1.0E-05; zero_muller ( func01, fatol, itmax, x1, x2, x3, xatol, xrtol, &xnew, &fxnew ); printf ( "\n" ); printf ( " X = %20.10f, %20.10f\n", creal ( xnew ), cimag ( xnew ) ); printf ( "\n" ); printf ( " with function value F(X):\n" ); printf ( "\n" ); printf ( " FX = %20.10f, %20.10f\n", creal ( fxnew ), cimag ( fxnew ) ); printf ( " ||FX|| = %20.10f\n", cabs ( fxnew ) ); return; } /******************************************************************************/ void test03 ( ) /******************************************************************************/ /* Purpose: test03() tests zero_muller() on Zhelyazkov's function Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt */ { double fatol; double complex fxnew; int itmax; int test; double complex x1; double complex x2; double complex x3; double complex xnew; double xatol; double xrtol; printf ( "\n" ); printf ( "test03():\n" ); printf ( " Demonstrate zero_muller() on Zhelyazkov''s function.\n" ); fatol = 1.0E-07; itmax = 10; for ( test = 1; test <= 2; test++ ) { /* First set of starting points. Result is X = ( 1.5705798926, 0.0 ) */ if ( test == 1 ) { x1 = 1.0 + 0.0 * I; x2 = 0.0 + 1.0 * I; x3 = 0.5 + 0.5 * I; } /* Second set of starting points. Result is X = ( -0.5802520567, 0.0 ). */ else if ( test == 2 ) { x1 = 0.0 + 1.0 * I; x2 = 1.0 + 2.0 * I; x3 = -1.0 + 2.0 * I; } xatol = 1.0E-07; xrtol = 1.0E-07; zero_muller ( func01, fatol, itmax, x1, x2, x3, xatol, xrtol, &xnew, &fxnew ); printf ( "\n" ); printf ( " X = %20.10f, %20.10f\n", creal ( xnew ), cimag ( xnew ) ); printf ( "\n" ); printf ( " with function value F(X):\n" ); printf ( "\n" ); printf ( " FX = %20.10f, %20.10f\n", creal ( fxnew ), cimag ( fxnew ) ); printf ( " ||FX|| = %20.10f\n", cabs ( fxnew ) ); } return; } /******************************************************************************/ double complex func01 ( double complex x ) /******************************************************************************/ /* Purpose: func01() evaluates F(X) = X*X+9. Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt Input: double complex X, the point at which the function is to be evaluated. Output: double complex FX, the function value at X. */ { double complex fx; fx = x * x + 9.0; return fx; } /******************************************************************************/ double complex func02 ( double complex x ) /******************************************************************************/ /* Purpose: func02() evaluates F(X) = (X*X+4)*(X-1)*(X+2). Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt Input: double complex X, the point at which the function is to be evaluated. Output: double complex FX, the function value at X. */ { double complex fx; fx = ( x * x + 4.0 ) * ( x - 10.0 ) * ( x + 20.0 ); return fx; } /******************************************************************************/ double complex func03 ( double complex z ) /******************************************************************************/ /* Purpose: func03() evaluates Zhelyazkov's function. Licensing: This code is distributed under the MIT license. Modified: 28 March 2024 Author: John Burkardt Input: double complex Z, the point at which the function is to be evaluated. Output: double complex FZ, the function value at Z. */ { double complex a; double complex b; double complex eps = { 0.4 + 0.0 * I }; double complex eta = { 0.64 + 0.0 * I }; double complex fz; double me = 0.384; double mo = 0.5; double complex of; double complex ok; double complex one = { 1.0 + 0.0 * I }; double x = 0.5; ok = z - me / csqrt ( eta ); of = z - mo; a = of * of + ( ok * ok ) * eta * ctanh ( x ); b = ( of - ok * eta ) / ( of - ok * eta * eta ); fz = of * of - one + ( eta * ok * ok - one ) * ctanh ( x ) - x * x * eps * eps * a * b; return fz; } /******************************************************************************/ void timestamp ( ) /******************************************************************************/ /* Purpose: timestamp() prints the current YMDHMS date as a time stamp. Example: 31 May 2001 09:45:54 AM Licensing: This code is distributed under the MIT license. Modified: 24 September 2003 Author: John Burkardt */ { # define TIME_SIZE 40 static char time_buffer[TIME_SIZE]; const struct tm *tm; time_t now; now = time ( NULL ); tm = localtime ( &now ); strftime ( time_buffer, TIME_SIZE, "%d %B %Y %I:%M:%S %p", tm ); fprintf ( stdout, "%s\n", time_buffer ); return; # undef TIME_SIZE }