Example 3: Jacobi

This example show the use of the "pcoord" intrinsic and the dot notation for acessing nearest neighbors.
/*
   This program implements the Jacobi algorithm to find the
   steady-state temperature distribution on an insulated
   two-dimensional plate, given constant boundary conditions.

   For SIZE = 128 and EPSILON = 0.1, the program should print:
      There are 11841 cells cooler than 50.00 degrees.

   For SIZE = 128 and EPSILON = 2.0, the program should print:
      There are 3852 cells cooler than 50.00 degrees

   For SIZE = 32 and EPSILON = 2.0, the program should print:
      There are 16044 cells cooler than 50.00 degrees (32 x 32 grid)
*/


#include <dpce.h>
#include <math.h>
#include <assert.h>
#include <stdio.h>

#define SIZE       32         /* Resolution of grid */

#define SHAPE_SIZE SIZE      /* Size of defined shape */

#define TEMP    50.0      /* Arbitrary cut-off  */
#define EPSILON 2.0

shape [SHAPE_SIZE][SHAPE_SIZE]cell;

int:cell  active;       /* 0 if cell boundary; 1 otherwise */
float:cell change;       /* Change in temperature */
float:cell new;          /* Newly calculated temperature */
float:cell old;          /* Previous temperature */



main () {
  float maxerr;          /* Largest change in temp. over grid */
  int   cool_cells;      /* Number of cells with temp. < TEMP */
  int   pass = 0;

     where ((pcoord(cell,0) < SIZE) && (pcoord(cell,1) < SIZE)) {
        where ((!pcoord(cell,0)) ||  (!pcoord(cell,1)) || 
               (pcoord(cell,1) == (SIZE-1))) 
        {
           old = new = 0.0;
           printf("A: %d positions active\n", += (int:cell) 1);
        } 
        else where (pcoord(cell,0) == (SIZE-1)) 
        {
           active = 0;
           old = new = 100.0;
           printf("B: %d positions active\n", += (int:cell) 1);
        } 
        else 
        {
           active = 1;
           new = 50.0;
           printf("C: %d positions active\n", += (int:cell) 1);
        }

        do {               /* Compute steady-state temperatures */
           where (active) {

	      /* BARR */
              old = new;

	      /* BARR */
              new = ([.-1][.]old + [.+1][.]old +
                     [.][.+1]old + [.][.-1]old) / 4.0;

              maxerr = >?= fabs(old-new);

              printf("Pass %3d, maxerr is %f\n", ++pass, maxerr);
           }
        } while (maxerr > EPSILON);

        cool_cells = (+= (new < TEMP));
     }  /* where */

  printf ("There are %d cells cooler than %5.2f degrees (%d x %d grid)\n",
     cool_cells, TEMP, SIZE, SIZE);

   if (SIZE == 32 && EPSILON == 2.0 && TEMP == 50.0)
      assert (cool_cells == 692);

   else if (SIZE == 128 && EPSILON == 0.1 && TEMP == 50.0)
      assert (cool_cells == 11841);

   else if (SIZE == 128 && EPSILON == 2.0 && TEMP == 50.0)
      assert (cool_cells == 3852);

   else
      {
      printf ("UNSUPPORTED #define COMBINATION FOR VALIDATION!\n");
      }
}

Options:

©1995 Pacific-Sierra Research Corporation. All rights reserved.

Send comments or suggestions to dpce2 at crescentbaysoftware.com.