Changes between Initial Version and Version 1 of cypress/Programming/HeatMassTransfer/ex23.c


Ignore:
Timestamp:
May 14, 2015 3:23:55 PM (4 years ago)
Author:
cmaggio
Comment:

migrated content from ccs wiki

Legend:

Unmodified
Added
Removed
Modified
  • cypress/Programming/HeatMassTransfer/ex23.c

    v1 v1  
     1= Ex23.c =
     2{{{#!C
     3/*
     4  Large Scale Computing
     5  Convective Heat/Mass Transfer
     6  ex23.c : Gauss-Seidel
     7 
     8  duphi/dx = d2phi/dx2
     9 */
     10#include<stdio.h>
     11#include<stdlib.h>
     12#include<math.h>
     13 
     14/* min and max macros */
     15#define max(a,b)        (((a) > (b)) ? (a) : (b))
     16#define min(a,b)        (((a) < (b)) ? (a) : (b))
     17 
     18int main(int argc, char **argv){
     19  double a,b,d,u;
     20  double *phi;
     21  double w;
     22  double dx,x,sol;
     23  double err;
     24  double ae,aw,ap;
     25  double tol = 1.0e-10;
     26  int i,num,itc;
     27  FILE *fp;
     28 
     29  if (argc < 6){
     30    printf("Usage:%s [NUM] [A] [B] [D] [U]\n",argv[0]);
     31    exit(-1);
     32  }
     33 
     34  /* set parameters */
     35  num = atoi(argv[1]);
     36  a = atof(argv[2]);
     37  b = atof(argv[3]);
     38  d = atof(argv[4]);
     39  u = atof(argv[5]);
     40 
     41  printf("num=%d A=%e B=%e D=%e U=%e\n",num,a,b,d,u);
     42 
     43  /* Memory Allocation and Check*/
     44  phi = (double *)calloc(num, sizeof(double));
     45  if (phi == (double *)NULL){
     46    printf("Memory Allocation Failed\n");
     47    exit(-1);   
     48  }
     49 
     50  /** Setup **/
     51  dx = 1.0 / (double)(num - 1);
     52  for (i = 0 ; i < num ; i++){
     53    phi[i] = 0.0;
     54  }
     55  ae = d / dx + max(-u,0.0);
     56  aw = d / dx + max(u,0.0);
     57  ap = ae + aw;
     58 
     59  /* Boundary Condition*/
     60  phi[0] = a;
     61  phi[num-1] = b;
     62 
     63  /* Solve with Gauss-Seidel Method*/
     64  itc = 0;
     65  while(1){
     66    /* update phi   &  Check Convergence */
     67    err = 0.0;   
     68    for (i = 1 ; i < (num - 1) ; i++){
     69      w =  (ae *  phi[i + 1] + aw *  phi[i - 1]) / ap;
     70      err += (w - phi[i]) * (w - phi[i]);     
     71      phi[i] = w;
     72    }
     73 
     74    itc++;
     75    if (err < tol) break;
     76    if ((itc % 100) == 0) printf("itc=%d err=%e\n",itc,err);
     77  }
     78 
     79  printf("Number of Iteration=%d\n",itc);
     80 
     81  /* Output Result */
     82  fp = fopen("res.dat","w");
     83  if (fp == NULL){
     84    printf("File could not create\n");
     85    exit(-1);
     86  }
     87 
     88  for (i = 0 ; i < num ; i++){
     89    x = dx * (double)i;
     90    sol = a + (b - a) * (exp(x * u / d) - 1.0) / (exp(u / d) - 1.0);
     91    fprintf(fp,"%e %e %e\n",x,phi[i],sol);
     92  }
     93  fclose(fp);
     94 
     95  /* END : Deallocate Memory */
     96  free(phi);
     97  printf("Done\n");
     98}
     99}}}