lukasz.proszek.info

Powrót do listy plików

---

backsubstitution.c

#include<stdio.h>
#include<math.h>
#include<ctype.h>
#define NWIER 50

typedef double typ;

int main(){
         double **U;
         double *b;
         double *x;
         double temp;
         int i,j;

        U = (typ **)calloc(NWIER,sizeof(typ));
        for(i=0; i<NWIER; i++){
                 U[i] = (typ *)calloc(NWIER,sizeof(typ));
        }//for i
   
        b= (typ **)calloc(NWIER,sizeof(typ));
        x= (typ **)calloc(NWIER,sizeof(typ));

        for(i=0;i<NWIER;i++){
                 for(j=0;j<NWIER;j++){
                                if(j>=i){
                                         U[i][j]=cos(.1*(i+j));
//                                       printf("U[%i][%i]=%f\n",i,j,U[i][j]);
                                }//if j>=i
                                else{
                                         U[i][j]=0;
                                }
                                
                 }//for j
        }//for i

        for(i=0;i<NWIER;i++){
                 x[i]=0;
                 b[i]=1;
        }
        
//      printf("%f,%f,%f",x[NWIER-1],b[NWIER-1],U[NWIER-1][NWIER-1]);
//      x[NWIER-1]=b[NWIER-1]/U[NWIER-1][NWIER-1];
//      ^- niepotrzebne bo x[49]=(b[49] - U[49][49]*(x[49]=0))/U[49][49]
//      czyli x[49]= b[49]/U[49][49]
        
//      printf("%f,%f,%f",x[NWIER-1],b[NWIER-1],U[NWIER-1][NWIER-1]);
//fortuna-mecukow 5.3.2 str 200
        for(j=NWIER-1;j>=0;j--){
                 temp=0.;
                 for(i=NWIER-1;i>j;i--){
                                temp+=U[j][i]*x[i];
                 }
                 x[j]=(b[j]-temp)/U[j][j];
        }
        
for(i=0;i<NWIER;i++)
         printf("%.10f\n",x[i]);

return 0;
}

---

syntax highlighted by Code2HTML, v. 0.9.1
kod przerobiony z HTML na XHTML za pomocą HTML Tidy for Linux/x86