/* University of Toronto -- CSC 270S/B70S -- Monday 11 January 1999 */ /* */ /* Homework Assignment 1: Root Finding */ /* Due: Friday 29 January at 4:00 pm -- Worth; 10% */ /******************************************************************** ** ** --- PLEASE FILL IN THIS SECTION COMPLETELY --- ** ** Last name: ** First name: ** student number: ** Username: ** Tutorial section: ** ********************************************************************/ /* * FILE: minmax.c * * Do not change the first part of this file. Look further down * to get to the part where you have to make changes. */ /******************************************************************** ** ** ** --- DO NOT MAKE CHANGES STARTING AT THIS LINE... ** ** **/ #include #include #include #include /* -- Symbolic constants. -- */ #define NUMBER_OF_FUNCTIONS 3 /* number of functions available */ #define PRINT_PRECISION 13 /* number of digits to print for */ /* floating-point values */ /* -- Type definitions. -- */ /* The type of value returned by the root-finding methods: */ /* needed so we can return two values with different types. */ typedef struct { long num_iter; /* number of iterations performed */ double root; /* approximate root found (if any) */ } result_type; /* -- Function prototypes (declarations). -- */ void print_usage( const char * ); int evaluate( int, double ); int find_extremum( int, char, double, double, double, long ); double f( int, double ); double df( int, double ); double ddf( int, double ); result_type bisection( int, double, double, double, long ); result_type fast_bisection( int, double, double, double, long ); result_type newton( int, double, double, long ); /* -- Main function. -- */ /* * FUNCTION: main * * This program can be called in one of two ways: by giving it a * function number and a real value (in which case it evaluates the * function and its derivatives at the given point) or by giving it a * function number, a method, an interval, a tolerance parameter, and * a maximum number of iterations (in which case it tries to find a * minimum or maximum of the function inside the given interval using * the specified method for root-finding). * * This main function simply checks the number of command-line * arguments and either calls the `evaluate' function or the * `find_extremum' function with the correctly parsed command-line * arguments. * * PARAMETERS: * int argc number of command-line arguments * char *argv[] array of command-line arguments * * RETURN VALUE: * int exit status (0 = OK) */ int main( int argc, char *argv[] ) { /* Based on the number of command-line arguments, call the */ /* appropriate function with the parsed arguments. */ if ( argc-1 == 2 ) { return evaluate( strtol(argv[1],NULL,0), strtod(argv[2],NULL) ); } else if ( argc-1 == 6 ) { return find_extremum( strtol(argv[1],NULL,0), tolower(argv[2][0]), strtod(argv[3],NULL), strtod(argv[4],NULL), strtod(argv[5],NULL), strtol(argv[6],NULL,0) ); } else { print_usage("Error: wrong number of arguments."); return argc; } /*-if-*/ } /*-main-*/ /* -- Function definitions. -- */ /* * FUNCTION: evaluate * * Evaluate the specified function and its derivatives at the given * point. * * ARGUMENTS: * int func_num the number of the function to evaluate * double value the real number at which to evaluate * * Preconditions: * NONE: The function is responsible for checking that the values * of the parameters are valid. * * RETURN VALUE: * int exit status (0 = OK, <0 = invalid parameter) * * Postconditions: * If the parameters are valid, the values of the specified * function and of its derivatives at the specified point are * printed on the screen. */ int evaluate( int func_num, double value ) { /* Check that the arguments are valid. */ if ( func_num < 1 || NUMBER_OF_FUNCTIONS < func_num ) { print_usage("Error: invalid function number."); return -1; } /*-if-*/ /* Output the values of the function and its derivatives. */ printf(" f%d(%.*g) = % .*g\n", func_num, PRINT_PRECISION, value, PRINT_PRECISION, f(func_num,value)); printf(" f%d'(%.*g) = % .*g\n", func_num, PRINT_PRECISION, value, PRINT_PRECISION, df(func_num,value)); printf(" f%d''(%.*g) = % .*g\n", func_num, PRINT_PRECISION, value, PRINT_PRECISION, ddf(func_num,value)); return 0; } /*-evaluate-*/ /* * FUNCTION: find_extremum * * Try to find a minimum or maximum of the specified function, using * the specified root-finding method to find a root of the function's * derivative. Output the results of the computation. * * ARGUMENTS: * int func_num the number of the function * char method the root-finding method to use * double x_left, x_right the starting interval * double epsilon the tolerance parameter * long max_iter the maximum number of iterations * * Preconditions: * NONE: The function is responsible for checking that the values * of the parameters are valid. * * RETURN VALUE: * int exit status (0 = OK, <0 = invalid parameter) * * Postconditions: * If the parameters are valid, the specified root-finding method * is called to find a zero of the specified function's derivative * inside the specified interval. The results of the computation * are printed on the screen. */ int find_extremum( int func_num, char method, double x_left, double x_right, double epsilon, long max_iter ) { result_type result; /* to store result of root-finding method */ /* Check that the arguments are valid. */ if ( func_num < 1 || NUMBER_OF_FUNCTIONS < func_num ) { print_usage("Error: invalid function number."); return -1; } /*-if-*/ if ( method != 'b' && method != 'f' && method != 'n' ) { print_usage("Error: invalid root-finding method."); return -1; } /*-if-*/ if ( x_left > x_right ) { print_usage("Error: invalid starting interval."); return -1; } /*-if-*/ if ( epsilon < 0 ) { print_usage("Error: invalid tolerance parameter."); return -1; } /*-if-*/ if ( max_iter < 0 ) { print_usage("Error: invalid maximum number of iterations."); return -1; } /*-if-*/ /* Output the values of the parameters, */ printf(" Finding an extremum of function number %d\n", func_num); printf(" in the interval [%.*g,%.*g]\n", PRINT_PRECISION, x_left, PRINT_PRECISION, x_right); printf(" with tolerance parameter %.*g\n", PRINT_PRECISION, epsilon); printf(" and maximum number of iterations %ld\n", max_iter); /* find a minimum or maximum (find a root of the derivative), */ switch ( method ) { case 'b': printf(" using the bisection method...\n"); result = bisection( func_num, x_left, x_right, epsilon, max_iter ); break; case 'f': printf(" using the \"fast bisection\" method...\n"); result = fast_bisection( func_num, x_left, x_right, epsilon, max_iter ); break; case 'n': printf(" using the Newton method...\n"); result = newton( func_num, ( x_left + x_right ) / 2.0, epsilon, max_iter ); break; default: } /*-switch-*/ /* and output the results. */ if ( result.num_iter < max_iter ) { printf("...found "); if ( ddf(func_num,result.root) < 0.0 ) { printf("a maximum"); } else if ( ddf(func_num,result.root) > 0.0 ) { printf("a minimum"); } else { printf("an extremum"); } /*-if-*/ printf(" v in %ld iterations:\n", result.num_iter); printf(" v = % .*g\n", PRINT_PRECISION, result.root); printf(" f(v) = % .*g\n", PRINT_PRECISION, f(func_num,result.root)); printf(" f'(v) = % .*g\n", PRINT_PRECISION, df(func_num,result.root)); printf(" f''(v) = % .*g\n", PRINT_PRECISION, ddf(func_num,result.root)); } else { printf("...found no extremum in %ld iterations.\n", result.num_iter); } /*-if-*/ return 0; } /*-find_extremum-*/ /* * FUNCTION: print_usage * * Print a description of how to invoke the program. * * PARAMETERS: * const char *error_message error message to print first * * Preconditions: * error_message != NULL * * RETURN VALUE: * void (NONE) * * Postconditions: * `error_message' is printed on the screen, followed by a * description of how to use the program. */ void print_usage( const char *error_message ) { puts(""); puts(error_message); puts(""); puts("To use the program, you must supply either two or six"); puts("command-line arguments. If you supply two arguments, they"); puts("must be a function number and a real value, in which case"); puts("the program will output the value of the given function and"); puts("of its derivatives at the given point. If you supply six"); puts("arguments, they must follow the following pattern:"); puts(""); puts(" minmax func_num method x_left x_right epsilon max_iter"); puts(""); puts("where:"); puts(""); puts(" - `func_num' is the number of the function for which you"); puts(" want to find a minimum or maximum (must be between 1 and"); printf(" %d);\n", NUMBER_OF_FUNCTIONS); puts(""); puts(" - `method' is the first letter of the root-finding method"); puts(" you want to use: 'b' for bisection, 'f' for \"fast"); puts(" bisection\", or 'n' for Newton;"); puts(""); puts(" - [x_left,x_right] is the starting interval inside which"); puts(" you want to find a minimum or maximum (x_left must be no"); puts(" greater than x_right);"); puts(""); puts(" - `epsilon' is the tolerance parameter, i.e., how small the"); puts(" interval must be before the root-finding method stops"); puts(" (must be nonnegative);"); puts(""); puts(" - `max_iter' is the maximum number of iterations to perform"); puts(" in the root-finding method (must be nonnegative)."); puts(""); puts("Try again!"); puts(""); } /*-print_usage-*/ /** ** ** ...DO NOT MAKE CHANGES STOPPING AT THIS LINE --- ** ** ** ********************************************************************/ /******************************************************************** ** ** ** Fill in the body of each of the functions below to perform ** ** the computation indicated in the comments. ** ** ** ********************************************************************/ /* * FUNCTION: f * * The function for which we want to find a minimum/maximum. * * PARAMETERS: * int func_num the number of the function to compute * double x the input to the function * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * * RETURN VALUE: * double the value of the function evaluated at x * * Postconditions: * The value of the function evaluated at x is returned. */ double f( int func_num, double x ) { switch ( func_num ) { case 1: return ( exp(x) - 3.5 * x * x - x + 5.001 ); case 2: return ( log(x) + sin( x * sqrt(x) ) ); case 3: /* Fill in code to compute the value of your function f. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computation. */ default: return 0.0; } /*-switch-*/ } /*-f-*/ /* * FUNCTION: df * * The derivative of the function above. * * PARAMETERS: * int func_num the number of the derivative to compute * double x the input to the derivative * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * * RETURN VALUE: * double the value of the derivative evaluated at x * * Postconditions: * The value of the derivative evaluated at x is returned. */ double df( int func_num, double x ) { switch ( func_num ) { case 1: return ( exp(x) - 7.0 * x - 1.0 ); case 2: return ( 1.0 / x + 1.5 * sqrt(x) * cos( x * sqrt(x) ) ); case 3: /* Fill in code to compute the derivative of your function f. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computation. */ default: return 0.0; } /*-switch-*/ } /*-df-*/ /* * FUNCTION: ddf * * The second derivative of the function above. * * PARAMETERS: * int func_num the number of the second derivative to compute * double x the input to the second derivative * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * * RETURN VALUE: * double the value of the second derivative evaluated at x * * Postconditions: * The value of the second derivative evaluated at x is returned. */ double ddf( int func_num, double x ) { switch ( func_num ) { case 1: return ( exp(x) - 7.0 ); case 2: return ( -1.0 / ( x * x ) + ( 0.75 / sqrt(x) ) * cos( x * sqrt(x) ) - 2.25 * x * sin( x * sqrt(x) ) ); case 3: /* Fill in code to compute the second derivative */ /* of your function f. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computation. */ default: return 0.0; } /*-switch-*/ } /*-ddf-*/ /* * FUNCTION: bisection * * The bisection method of root-finding: finds an approximate * root of the derivative of the given function. * * PARAMETERS: * int func_num the number of the function f * double x_left, x_right the interval inside which the root lies * double epsilon the tolerance parameter * long max_iter the maximum number of iterations * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * x_left <= x_right * df(func_num,x_left) * df(func_num,x_right) < 0.0 * df is continuous on the interval [x_left,x_right] * 0.0 <= epsilon * 0 <= max_iter * * RETURN VALUE: * result_type * .num_iter the number of iterations performed * .root an approximate root of df inside [x_left,x_right] * * Postconditions: * the function returns within at most `max_iter' iterations * if the function returns before `max_iter' iterations have been * completed, then the approximate root is within `epsilon' of a * real root */ result_type bisection( int func_num, double x_left, double x_right, double epsilon, long max_iter ) { result_type result; long num_iter = 0; double root = 0.0; /* Fill in the rest of the `bisection' function. Use */ /* `num_iter' to keep track of the number of iterations, and */ /* make sure to set `root' equal to the approximate root found. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computations. */ /* This stores the number of iterations and the approximate */ /* root in a structure that can be returned by the function. */ result.num_iter = num_iter; result.root = root; return result; } /*-bisection-*/ /* * FUNCTION: fast_bisection * * The "fast bisection" method of root-finding: finds an approximate * root of the derivative of the given function. * * PARAMETERS: * int func_num the number of the function f * double x_left, x_right the interval inside which the root lies * double epsilon the tolerance parameter * long max_iter the maximum number of iterations * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * x_left <= x_right * df(func_num,x_left) * df(func_num,x_right) < 0.0 * df is continuous on the interval [x_left,x_right] * 0.0 <= epsilon * 0 <= max_iter * * RETURN VALUE: * result_type * .num_iter the number of iterations performed * .root an approximate root of df inside [x_left,x_right] * * Postconditions: * the function returns within at most `max_iter' iterations * if the function returns before `max_iter' iterations have been * completed, then the approximate root is within `epsilon' of a * real root */ result_type fast_bisection( int func_num, double x_left, double x_right, double epsilon, long max_iter ) { result_type result; long num_iter = 0; double root = 0.0; /* Fill in the rest of the `fast_bisection' function. Use */ /* `num_iter' to keep track of the number of iterations, and */ /* make sure to set `root' equal to the approximate root found. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computations. */ /* This stores the number of iterations and the approximate */ /* root in a structure that can be returned by the function. */ result.num_iter = num_iter; result.root = root; return result; } /*-fast_bisection-*/ /* * FUNCTION: newton * * The Newton method of root-finding: finds an approximate * root of the derivative of the given function. * * PARAMETERS: * int func_num the number of the function f * double x0 the initial guess for a root * double epsilon the tolerance parameter * long max_iter the maximum number of iterations * * Preconditions: * 1 <= func_num <= NUMBER_OF_FUNCTIONS * df is continuous on some interval around x0 * ddf is continuous and non-zero on that same interval * 0.0 <= epsilon * 0 <= max_iter * * RETURN VALUE: * result_type * .num_iter the number of iterations performed * .root an approximate root of df * * Postconditions: * the function returns within at most `max_iter' iterations * if the function returns before `max_iter' iterations have been * completed, then the approximate root is within `epsilon' of a * real root */ result_type newton( int func_num, double x0, double epsilon, long max_iter ) { result_type result; long num_iter = 0; double root = 0.0; /* Fill in the rest of the `newton' function. Use */ /* `num_iter' to keep track of the number of iterations, and */ /* make sure to set `root' equal to the approximate root found. */ /* Do NOT change the code that is already here, only add what */ /* is needed to perform the necessary computations. */ /* This stores the number of iterations and the approximate */ /* root in a structure that can be returned by the function. */ result.num_iter = num_iter; result.root = root; return result; } /*-newton-*/