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

Timestamp:

05/14/15 15:19:00 (10 years ago)

Author:

cmaggio

Comment:

migrated content from ccs wiki

cypress/Programming/HeatMassTransfer/ex17.c

@@ +1 @@
+= Ex17.c =
+{{{#!C
+/*
+  Large Scale Computing
+  2D Heat/Mass Transfer
+  ex17.c : Numerical Analysis
+  Solve for
+  d2c/dx2 + d2c/dy2 + 1 = 0
+  With the boundary conditions of c = 0 
+  along lines of x=1,-1, and y=1, and -1
+ */
+#include<stdio.h>
+#include<stdlib.h>
+#include<math.h>
+ 
+#define NUM 20  /* number of points for x and y*/
+#define TOL 1.0e-12  /* convergence tolerance */
+ 
+int main(void ){
+  int i,j,k,n,ie,iw,in,is;
+  double x,y,dx,w,err;
+  double cnc[NUM * NUM]; /* define 1D array of length NUM * NUM */
+  FILE *fp;
+ 
+  /*assuming dx = dy : domain is 2x2 */
+  dx = 2.0 / (double)(NUM - 1);
+ 
+  /* Initialize */
+  for (n = 0 ; n < NUM * NUM ; n++){
+    cnc[n] = 0.0;
+  }
+ 
+  /* computing for C with Gauss-Seidel */
+  k = 0;
+  while(1){
+    err = 0.0;
+    for (n = 0 ; n < NUM * NUM ; n++){
+      i = n % NUM; /* i is residual */
+      j = n / NUM; /* j is division */
+      if (i == 0) continue; /*skip left end*/
+      if (i == NUM-1) continue; /*skip right end*/
+      if (j == 0) continue; /*skip bottom end*/
+      if (j == NUM-1) continue; /*skip top end*/
+ 
+      ie = n + 1;
+      iw = n - 1;
+      in = n + NUM;
+      is = n - NUM;
+      w = 0.25 * (dx * dx + cnc[ie] + cnc[iw] + cnc[in] + cnc[is]);
+      err += (w - cnc[n]) * (w - cnc[n]);
+      cnc[n] = w;
+    }
+    k++;
+    if (err < TOL) break;
+  }
+  printf("# of Iteration=%d\n",k);
+ 
+  /* Output Result */
+  fp = fopen("res.dat","w");
+  if (fp == NULL){
+    printf("File could not create\n");
+    exit(-1);
+  }
+ 
+  for (i = 0 ;i < NUM ; i++){
+    for (j = 0 ; j < NUM ; j++){
+      x = -1.0 + 2.0 * (double)i / (double)(NUM - 1);
+      y = -1.0 + 2.0 * (double)j / (double)(NUM - 1);
+      fprintf(fp,"%e %e %e\n",x,y,cnc[i + NUM * j]);
+    }
+  }
+  fclose(fp);
+ 
+  printf("Done\n");
+}
+}}}