C語言 例題4-6 # 利用 辛普森 理則計算 雙重積分 f(x,y)= x^2 y 在[1,2] dx 與 [0,1] dy 的定積分

C語言 例題4-6 # 利用 辛普森 理則計算 雙重積分 f(x,y)= x^2 y 在[1,2] dx 與 [0,1] dy 的定積分

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define PI 3.141596
#define N 2

// This function multiplies
// mat1[][] and mat2[][], and
// stores the result in res[][]
void multiply(double mat1[N+1],double mat2[N+1],double res[][N+1])
{
    int i, j;
    for (i = 0; i < N+1; i++)
    {
        for (j = 0; j < N+1; j++)
        {
            res[i][j] = 0;
            res[i][j] += mat1[i]*mat2[j];
        }
    }
}

void multiply2(double res[N+1][N+1],double ff[N+1][N+1],double ffm[N+1][N+1])
{
    int i, j;
    for (i = 0; i < N+1; i++)
    {
        for (j = 0; j < N+1; j++)
        {
            ffm[i][j] = 0;
            ffm[i][j] += res[i][j]*ff[i][j];
        }
    }
}


// Driver Code
int main()
{
    int i, j ,n=2;
    double res[N+1][N+1]; // To store result
    double mat1[N+1] = {1, 4, 1};
    double mat2[N+1] = {1, 4, 1};
    double x[N+1], y[N+1], ff[N+1][N+1] ,ffm[N+1][N+1];
 
    multiply(mat1, mat2, res);

    printf("  Result matrix is \n");
    for (i = 0; i <= N; i++)
    {
        for (j = 0; j <=N; j++)
            printf(" %0.2lf " ,res[i][j]);
        printf("\n\n");
    }
 
    double a=1 ,b=2, hx ;
    hx=(b-a)/n;
    x[0]=a;
    printf("x[] = %0.2lf ",x[0]);
    for (i = 1; i <= N; i++)
    {
        x[i]=i*hx+a;
        printf(" %0.2lf ",x[i]);
    }
    printf("\n");
 
    printf("x^2 =");
    for (i = 0; i <= N; i++)
    {
        x[i]*=x[i];           //X^2 = X*X 
        printf(" %0.2lf ",x[i]);
    }
    printf("\n");



    double c=0 ,d=1, hy ;
    hy=(d-c)/n;
    y[0]=c;
    printf("y[] = %0.2lf ",y[0]);
    for (i = 1; i <= N; i++)
    {
        y[i]=i*hy+c;
        printf(" %0.2lf ",y[i]);
    }
    printf("\n");

    multiply(x , y , ff);
 
    printf("\n\n  Result matrix is \n");
    for (i = 0; i <= N; i++)
    {
        for (j = 0; j <=N; j++)
            printf(" %0.4lf " ,ff[i][j]);
        printf("\n\n");
    }
    multiply2(res , ff , ffm);
 
    double sum=0.0;
    printf("\n\n  Result matrix is \n");
    for (i = 0; i <= N; i++)
    {
        for (j = 0; j <=N; j++)
            {
                sum=sum+ffm[i][j];
                printf(" %0.4lf " ,ffm[i][j]);
            } 
        printf("\n\n");
    }
    sum=(1.0/9)*hy*hx*sum;
    printf("雙重積分值= %0.8lf",sum);
 
    return 0;
}

輸出畫面
  Result matrix is
 1.00  4.00  1.00

 4.00  16.00  4.00

 1.00  4.00  1.00

x[] = 1.00  1.50  2.00
x^2 = 1.00  2.25  4.00
y[] = 0.00  0.50  1.00


  Result matrix is
 0.0000  0.5000  1.0000

 0.0000  1.1250  2.2500

 0.0000  2.0000  4.0000



  Result matrix is
 0.0000  2.0000  1.0000

 0.0000  18.0000  9.0000

 0.0000  8.0000  4.0000

雙重積分值= 1.16666667