Subversion Repositories programming

Rev

Rev 156 | Rev 158 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 156 Rev 157
Line 9... Line 9...
9 
 * - Implemented Functions 1-3 from the Project #3 Handout.
 9 
 * - Implemented Functions 1-3 from the Project #3 Handout.
10 
 * - Implemented Algorithm 7 (LaGrange) from Notes Set #3.
 10 
 * - Implemented Algorithm 7 (LaGrange) from Notes Set #3.
11 
 * - Implemented Algorithm 8-9 (Newton Divided Diff) from Notes Set #3.
 11 
 * - Implemented Algorithm 8-9 (Newton Divided Diff) from Notes Set #3.
12 
 * - Implemented GenEvenPts() which will generate evenly spaced points.
 12 
 * - Implemented GenEvenPts() which will generate evenly spaced points.
13 
 * - Implemented GenChebychevPts() which will generate Chebychev points.
 13 
 * - Implemented GenChebychevPts() which will generate Chebychev points.
- 
 
 14 
 * - Split the Newton algorithm into two parts: NewtonMethod() and NewtonCoef().
14 
 *
 15 
 *
15 
 */
 16 
 */
16 
 
 17 
 
17 
#include <cstdio>
 18 
#include <cstdio>
18 
#include <cmath>
 19 
#include <cmath>
19 
using namespace std;
 20 
using namespace std;
20 
 
 21 
 
21 
#define X_MIN = -5.0
 22 
#define MAX_POINTS 26
- 
 
 23 
 
- 
 
 24 
// Global Variables for easier return values
22 
#define X_MAX = 5.0
 25 
double absmax, xsave;
23 
 
 26 
 
24 
/**
 27 
/**
25 
 * Function #1 from the Project #3 Handout.
 28 
 * Function #1 from the Project #3 Handout.
26 
 *
 29 
 *
27 
 * @param x the place to calculate the value of this function
 30 
 * @param x the place to calculate the value of this function
Line 83... Line 86...
83 
 
 86 
 
84 
    return value;
 87 
    return value;
85 
}
 88 
}
86 
 
 89 
 
87 
/**
 90 
/**
88 
 * Find the value of the Newton Interpolating Polynomial at a point.
 91 
 * Find the Newton Coefficients.
89 
 * The a[i]'s are found using a divided difference table.
 92 
 * This only needs to be called once per set of x,y values.
90 
 *
 93 
 *
91 
 * @param x the x[i] values
 94 
 * @param x the x[i] values
92 
 * @param y the y[i] values
 95 
 * @param y the y[i] values
- 
 
 96 
 * @param a the a[i] values will be returned here
93 
 * @param n the number of points
 97 
 * @param n the number of points
94 
 * @param point the point to evaluate at
 - 
 
95 
 * @return the value of the Interpolating Poly at point
 - 
 
96 
 */
 98 
 */
97 
double NewtonMethod (double *x, double *y, int n, double point)
 99 
void NewtonCoef (double *x, double *y, double *a, int n)
98 
{
 100 
{
99 
    int i, j;
 101 
    int i, j;
100 
    double a[n];
 - 
 
101 
 
 102 
 
102 
    for (i=0; i<n; i++)
 103 
    for (i=0; i<n; i++)
103 
        a[i] = y[i];
 104 
        a[i] = y[i];
104 
 
 105 
 
105 
    for (i=1; i<n; i++)
 106 
    for (i=1; i<n; i++)
106 
        for (j=n-1; j>=i; j--)
 107 
        for (j=n-1; j>=i; j--)
107 
            a[j] = (a[j] - a[j-1]) / (x[j] - x[j-i]);
 108 
            a[j] = (a[j] - a[j-1]) / (x[j] - x[j-i]);
- 
 
 109 
}
108 
 
 110 
 
- 
 
 111 
/**
109 
    // At this point, all of the a[i]'s have been calculated,
 112 
 * Find the value of the Newton Interpolating Polynomial at a point.
110 
    // using a Divided Difference Table.
 113 
 * The a[i]'s are found using a divided difference table.
- 
 
 114 
 *
- 
 
 115 
 * @param x the x[i] values
- 
 
 116 
 * @param y the y[i] values
- 
 
 117 
 * @param a the a[i] values (Newton Coefficients)
- 
 
 118 
 * @param n the number of points
- 
 
 119 
 * @param point the point to evaluate at
- 
 
 120 
 * @return the value of the Interpolating Poly at point
- 
 
 121 
 */
- 
 
 122 
double NewtonMethod (double *x, double *y, double *a, int n, double point)
111 
 
 123 
{
- 
 
 124 
    int i;
112 
    double xterm = 1.0;
 125 
    double xterm = 1.0;
113 
    double value = 0.0;
 126 
    double value = 0.0;
114 
 
 127 
 
115 
    for (i=0; i<n; i++)
 128 
    for (i=0; i<n; i++)
116 
    {
 129 
    {
Line 162... Line 175...
162 
        x[i] = -5.0 * cos((float)i * M_PI / (float)(num_pts-1));
 175 
        x[i] = -5.0 * cos((float)i * M_PI / (float)(num_pts-1));
163 
        y[i] = func(x[i]);
 176 
        y[i] = func(x[i]);
164 
    }
 177 
    }
165 
}
 178 
}
166 
 
 179 
 
- 
 
 180 
/**
- 
 
 181 
 * Brute-force search the Newton Interpolating Poly for the x and y values.
- 
 
 182 
 * The error calculations will be for the given function (func1, func2, or func3).
- 
 
 183 
 *
- 
 
 184 
 * NOTE: The "return values" are stored in two global variables: absmax and xsave.
- 
 
 185 
 *
- 
 
 186 
 * @param x the x[i] values
- 
 
 187 
 * @param y the y[i] values
- 
 
 188 
 * @param n the number of points for this polynomial
- 
 
 189 
 * @param func the function to use for the error calculation
- 
 
 190 
 */
- 
 
 191 
void brute_search_newton (double *x, double *y, int n, double(*func)(double))
- 
 
 192 
{
- 
 
 193 
    absmax = 0.0;
- 
 
 194 
    double xval = -5.0;
- 
 
 195 
    double h = 10.0/500.0;
- 
 
 196 
    double temp, error;
- 
 
 197 
    double a[n];
- 
 
 198 
    int i;
- 
 
 199 
 
- 
 
 200 
    NewtonCoef(x, y, a, n);
- 
 
 201 
 
- 
 
 202 
    for (i=0; i<=500; i++)
- 
 
 203 
    {
- 
 
 204 
        temp = NewtonMethod(x, y, a, n, xval);
- 
 
 205 
        error = abs(temp - func(xval));
- 
 
 206 
 
- 
 
 207 
        if (error > absmax)
- 
 
 208 
        {
- 
 
 209 
            absmax = error;
- 
 
 210 
            xsave = xval;
- 
 
 211 
        }
- 
 
 212 
 
- 
 
 213 
        xval += h;
- 
 
 214 
    }
- 
 
 215 
}
- 
 
 216 
 
- 
 
 217 
/**
- 
 
 218 
 * Brute-force search the LaGrange Interpolating Poly for the x and y values.
- 
 
 219 
 * The error calculations will be for the given function (func1, func2, or func3).
- 
 
 220 
 *
- 
 
 221 
 * NOTE: The "return values" are stored in two global variables: absmax and xsave.
- 
 
 222 
 *
- 
 
 223 
 * @param x the x[i] values
- 
 
 224 
 * @param y the y[i] values
- 
 
 225 
 * @param n the number of points for this polynomial
- 
 
 226 
 * @param func the function to use for the error calculation
- 
 
 227 
 */
- 
 
 228 
void brute_search_lagrange (double *x, double *y, int n, double(*func)(double))
- 
 
 229 
{
- 
 
 230 
    absmax = 0.0;
- 
 
 231 
    double xval = -5.0;
- 
 
 232 
    double h = 10.0/500.0;
- 
 
 233 
    double temp, error;
- 
 
 234 
    int i;
- 
 
 235 
 
- 
 
 236 
    for (i=0; i<=500; i++)
- 
 
 237 
    {
- 
 
 238 
        temp = LaGrangeMethod(x, y, n, xval);
- 
 
 239 
        error = abs(temp - func(xval));
- 
 
 240 
 
- 
 
 241 
        if (error > absmax)
- 
 
 242 
        {
- 
 
 243 
            absmax = error;
- 
 
 244 
            xsave = xval;
- 
 
 245 
        }
- 
 
 246 
 
- 
 
 247 
        xval += h;
- 
 
 248 
    }
- 
 
 249 
}
- 
 
 250 
 
167 
int main (void)
 251 
int main (void)
168 
{
 252 
{
169 
    const int size = 11;
 253 
    int size, i;
170 
    double x[size];
 254 
    double x[MAX_POINTS];
171 
    double y[size];
 255 
    double y[MAX_POINTS];
- 
 
 256 
 
- 
 
 257 
    // In the following code, the correct functions are run for each of the
- 
 
 258 
    // 12 times this needs to be run.
- 
 
 259
- 
 
 260 
    printf ("Newton Polynomial, Equal Points, Function #1\n");
- 
 
 261 
    printf ("============================================================\n");
- 
 
 262 
    for (i=6; i<=26; i+=5)
- 
 
 263 
    {
- 
 
 264 
        size = i;
- 
 
 265 
        GenEvenPts (size, &func1, x, y);
- 
 
 266 
        brute_search_newton (x, y, size, &func1);
- 
 
 267 
 
- 
 
 268 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 269 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 270 
    }
- 
 
 271
- 
 
 272 
    printf ("\n\n");
- 
 
 273 
    printf ("Newton Polynomial, Equal Points, Function #2\n");
- 
 
 274 
    printf ("============================================================\n");
- 
 
 275 
    for (i=6; i<=26; i+=5)
- 
 
 276 
    {
- 
 
 277 
        size = i;
- 
 
 278 
        GenEvenPts (size, &func2, x, y);
- 
 
 279 
        brute_search_newton (x, y, size, &func2);
- 
 
 280 
 
- 
 
 281 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 282 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 283 
    }
- 
 
 284
- 
 
 285 
    printf ("\n\n");
- 
 
 286 
    printf ("Newton Polynomial, Equal Points, Function #3\n");
- 
 
 287 
    printf ("============================================================\n");
- 
 
 288 
    for (i=6; i<=26; i+=5)
- 
 
 289 
    {
- 
 
 290 
        size = i;
- 
 
 291 
        GenEvenPts (size, &func3, x, y);
- 
 
 292 
        brute_search_newton (x, y, size, &func3);
- 
 
 293 
 
- 
 
 294 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 295 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 296 
    }
- 
 
 297
- 
 
 298 
    printf ("\n\n");
- 
 
 299 
    printf ("LaGrange Polynomial, Equal Points, Function #1\n");
- 
 
 300 
    printf ("============================================================\n");
- 
 
 301 
    for (i=6; i<=26; i+=5)
- 
 
 302 
    {
- 
 
 303 
        size = i;
- 
 
 304 
        GenEvenPts (size, &func1, x, y);
- 
 
 305 
        brute_search_lagrange (x, y, size, &func1);
- 
 
 306 
 
- 
 
 307 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 308 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 309 
    }
- 
 
 310
- 
 
 311 
    printf ("\n\n");
- 
 
 312 
    printf ("LaGrange Polynomial, Equal Points, Function #2\n");
- 
 
 313 
    printf ("============================================================\n");
- 
 
 314 
    for (i=6; i<=26; i+=5)
- 
 
 315 
    {
172 
    double point = 2.0;
 316 
        size = i;
- 
 
 317 
        GenEvenPts (size, &func2, x, y);
- 
 
 318 
        brute_search_lagrange (x, y, size, &func2);
173 
 
 319 
 
- 
 
 320 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 321 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 322 
    }
- 
 
 323
- 
 
 324 
    printf ("\n\n");
- 
 
 325 
    printf ("LaGrange Polynomial, Equal Points, Function #3\n");
- 
 
 326 
    printf ("============================================================\n");
- 
 
 327 
    for (i=6; i<=26; i+=5)
- 
 
 328 
    {
- 
 
 329 
        size = i;
174 
    //GenEvenPts (size, &func1, x, y);
 330 
        GenEvenPts (size, &func3, x, y);
175 
    GenChebychevPts (size, &func1, x, y);
 331 
        brute_search_lagrange (x, y, size, &func3);
176 
 
 332 
 
- 
 
 333 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 334 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 335 
    }
- 
 
 336
- 
 
 337 
    printf ("\n\n");
- 
 
 338 
    printf ("Newton Polynomial, Chebychev Points, Function #1\n");
- 
 
 339 
    printf ("============================================================\n");
177 
    for (int i=0; i<size; i++)
 340 
    for (i=6; i<=26; i+=5)
- 
 
 341 
    {
- 
 
 342 
        size = i;
- 
 
 343 
        GenChebychevPts (size, &func1, x, y);
178 
        printf("x[%d] = %e -- y[%d] = %e\n", i, x[i], i, y[i]);
 344 
        brute_search_newton (x, y, size, &func1);
179 
 
 345 
 
- 
 
 346 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 347 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 348 
    }
- 
 
 349
- 
 
 350 
    printf ("\n\n");
180 
    printf ("LaGrange = %e\n", LaGrangeMethod(x, y, size, point));
 351 
    printf ("Newton Polynomial, Chebychev Points, Function #2\n");
- 
 
 352 
    printf ("============================================================\n");
- 
 
 353 
    for (i=6; i<=26; i+=5)
- 
 
 354 
    {
- 
 
 355 
        size = i;
- 
 
 356 
        GenChebychevPts (size, &func2, x, y);
- 
 
 357 
        brute_search_newton (x, y, size, &func2);
- 
 
 358 
 
- 
 
 359 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 360 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 361 
    }
- 
 
 362
- 
 
 363 
    printf ("\n\n");
181 
    printf ("Newton = %e\n", NewtonMethod(x, y, size, point));
 364 
    printf ("Newton Polynomial, Chebychev Points, Function #3\n");
- 
 
 365 
    printf ("============================================================\n");
- 
 
 366 
    for (i=6; i<=26; i+=5)
- 
 
 367 
    {
- 
 
 368 
        size = i;
- 
 
 369 
        GenChebychevPts (size, &func3, x, y);
- 
 
 370 
        brute_search_newton (x, y, size, &func3);
182 
 
 371 
 
- 
 
 372 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 373 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 374 
    }
- 
 
 375
- 
 
 376 
    printf ("\n\n");
- 
 
 377 
    printf ("LaGrange Polynomial, Chebychev Points, Function #1\n");
- 
 
 378 
    printf ("============================================================\n");
- 
 
 379 
    for (i=6; i<=26; i+=5)
- 
 
 380 
    {
- 
 
 381 
        size = i;
- 
 
 382 
        GenChebychevPts (size, &func1, x, y);
- 
 
 383 
        brute_search_lagrange (x, y, size, &func1);
- 
 
 384 
 
- 
 
 385 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 386 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 387 
    }
- 
 
 388
- 
 
 389 
    printf ("\n\n");
- 
 
 390 
    printf ("LaGrange Polynomial, Chebychev Points, Function #2\n");
- 
 
 391 
    printf ("============================================================\n");
- 
 
 392 
    for (i=6; i<=26; i+=5)
- 
 
 393 
    {
- 
 
 394 
        size = i;
- 
 
 395 
        GenChebychevPts (size, &func2, x, y);
- 
 
 396 
        brute_search_lagrange (x, y, size, &func2);
- 
 
 397 
 
- 
 
 398 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 399 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 400 
    }
- 
 
 401
- 
 
 402 
    printf ("\n\n");
- 
 
 403 
    printf ("LaGrange Polynomial, Chebychev Points, Function #3\n");
- 
 
 404 
    printf ("============================================================\n");
- 
 
 405 
    for (i=6; i<=26; i+=5)
- 
 
 406 
    {
- 
 
 407 
        size = i;
- 
 
 408 
        GenChebychevPts (size, &func3, x, y);
- 
 
 409 
        brute_search_lagrange (x, y, size, &func3);
- 
 
 410 
 
- 
 
 411 
        printf ("For %-2d Points, MAX ERROR is: %e"
- 
 
 412 
                " -- OCCURS at x = %e\n", size, absmax, xsave);
- 
 
 413 
    }
- 
 
 414
183 
    return 0;
 415 
    return 0;
184 
}
 416 
}
185 
 
 417