Changeset 459 for trunk/matml

Timestamp:

06/17/2009 03:52:54 PM (3 years ago)

Author:

powell

Message:

Major update splitting out two libraries and updating documentation.

Location:

trunk/matml/src/ternary

Files:

• book.c (modified) (9 diffs)
• ChangeLog (modified) (2 diffs)
• configure.in (modified) (1 diff)
• freenergy.c (modified) (3 diffs)
• geomview.c (modified) (6 diffs)
• Makefile.am (modified) (3 diffs)
• qhull.c (modified) (14 diffs)
• spinodal.c (modified) (10 diffs)
• square.c (modified) (7 diffs)
• ternary.c (modified) (6 diffs)
• ternary.h (modified) (8 diffs)
• ternary.tex.in (modified) (2 diffs)

trunk/matml/src/ternary/book.c

@@ -8 +8 @@
 
 
+#include "freenergy.h" /*+ Free energy prototypes, typedefs, etc. +*/
 #include "ternary.h" /*+ Ternary prototypes, typedefs, etc. +*/
 #include <math.h>    /*+ For sqrt(), fmin()/fmax() +*/
+#include <stdlib.h>  /*+ For malloc and free +*/
 
 
@@ -32 +34 @@
   side).
 
-  ternary_point *points Array of ternary points of size
+  energy_point *points Array of energy points of size
   +latex+$(N+1)(N+2)/2$, where $N$ is the resolution.
   -latex-(N+1)(N+2)/2, where N is the resolution.
@@ -53 +55 @@
 
 int init_triangle_array
-(int resolution, ternary_point *points, int efunc,
+(int resolution, energy_point *points, int efunc,
  double A2, double A3, double B2, double B3, double C2, double C3)
 {
@@ -107 +109 @@
   -latex-C3-direction.
 
-  ternary_point *points Array of ternary points of size
+  energy_point *points Array of energy points of size
   +latex+({\tt res2}+1)$\times$({\tt res3}+1).
   -latex-(res2+1) * (res3+1).
@@ -132 +134 @@
 
 int init_rectangle_array
-(int res2, int res3, ternary_point *points, int efunc,
+(int res2, int res3, energy_point *points, int efunc,
  double A2, double A3, double B2, double B3)
 {
@@ -243 +245 @@
   and colors the points according to free energy as well.
 
-  int scale_energy_array It returns zero or an error code.
+  int energy_to_visual_array It returns zero or an error code.
 
   int num_points Number of points in the array
 
-  ternary_point *points Array of ternary points of size
+  energy_point *epoints Array of energy points of size
   +latex+$(N+1)(N+2)/2$, where $N$ is the resolution.
   -latex-(N+1)(N+2)/2, where N is the resolution.
+
+  visual_point **tpoints Pointer to array of ternary points of size
+  +latex+$(N+1)(N+2)/2$, where $N$ is the resolution.
+  -latex-(N+1)(N+2)/2, where N is the resolution.
+  This will free and re-create this array if present.
 
   double y0 Reference (minimum)
@@ -277 +284 @@
   ++++++++++++++++++++++++++++++++++++++*/
 
-int scale_energy_array (int num_points, ternary_point *points,
-                        double y0, double G0, double y1, double G1,
-                        double *Gmin, double *Gmax, int square)
+int energy_to_visual_array
+(int num_points, energy_point *epoints, visual_point **tpoints, double y0,
+ double G0, double y1, double G1, double *Gmin, double *Gmax, int square)
 {
   int i;
 
+  if (*tpoints)
+    free (*tpoints);
+
+  if (!(*tpoints = (visual_point *) malloc (num_points * sizeof (visual_point))))
+    return -1;
+
   if (G0==0. && G1==0.)
     {
-      G0 = G1 = points[0].G;
+      G0 = G1 = epoints[0].G;
       for (i=0; i<num_points; i++)
         {
-          G0 = fmin (points[i].G, G0);
-          G1 = fmax (points[i].G, G1);
+          G0 = fmin (epoints[i].G, G0);
+          G1 = fmax (epoints[i].G, G1);
         }
       if (Gmin)
@@ -300 +313 @@
       if (Gmin)
         {
-          *Gmin = points[0].G;
+          *Gmin = epoints[0].G;
           for (i=0; i<num_points; i++)
-            *Gmin = fmin (points[i].G, *Gmin);
+            *Gmin = fmin (epoints[i].G, *Gmin);
         }
       if (Gmax)
         {
-          *Gmax = points[0].G;
+          *Gmax = epoints[0].G;
           for (i=0; i<num_points; i++)
-            *Gmax = fmax (points[i].G, *Gmax);
+            *Gmax = fmax (epoints[i].G, *Gmax);
         }
     }
@@ -321 +334 @@
       if (square)
         {
-          points[i].x = points[i].C2;
-          points[i].z = 1.-points[i].C3;
+          (*tpoints)[i].x = epoints[i].C2;
+          (*tpoints)[i].z = 1.-epoints[i].C3;
         }
       else
         {
-          points[i].x = points[i].C2 + 0.5*points[i].C3;
-          points[i].z = (1. - points[i].C3) * sqrt(3)/2.;
-        }
-      points[i].y = Grel(points[i].G);
-      points[i].red   = RED   (Grel (points[i].G));
-      points[i].green = GREEN (Grel (points[i].G));
-      points[i].blue  = BLUE  (Grel (points[i].G));
-      points[i].alpha = 1.;
-    }
-
-  return 0;
-}
+          (*tpoints)[i].x = epoints[i].C2 + 0.5*epoints[i].C3;
+          (*tpoints)[i].z = (1. - epoints[i].C3) * sqrt(3)/2.;
+        }
+      (*tpoints)[i].y = Grel(epoints[i].G);
+      (*tpoints)[i].red   = RED   (Grel (epoints[i].G));
+      (*tpoints)[i].green = GREEN (Grel (epoints[i].G));
+      (*tpoints)[i].blue  = BLUE  (Grel (epoints[i].G));
+      (*tpoints)[i].alpha = 1.;
+    }
+
+  return 0;
+}

trunk/matml/src/ternary/ChangeLog

@@ -3 +3 @@
   * Overhauled the free energy interface and functions, adding derivative
     implementations.
+  * Added WITH_DERIVATIVES and NO_INFINITY options to free energies().
+  * Split free energy out into its own library, for use in phase field codes.
   * Added phase_boundary object.
   * Changed hull calculation interface to only hold the qhull variables during
@@ -9 +11 @@
     minimization algorithm.
   * Report all phase boundaries in the hull as phase_boundary objects.
-  * Added WITH_DERIVATIVES and NO_INFINITY options to free energies().
   * Added spinodal.c with calc_spinodal function returning a phase_boundary.
   * New "square" free energy parameters for pseudo-square representation.
   * Rectangle bookkeeping functions complete square functionality.
+  * Split all non-free-energy functions out into a ternary library.
   * New "square" example program demonstrates these functions.
 

trunk/matml/src/ternary/configure.in

@@ -18 +18 @@
 AC_C_INLINE
 AM_PROG_CC_C_O
+AM_PROG_LIBTOOL
 
 dnl Geomview check

trunk/matml/src/ternary/freenergy.c

@@ -6 +6 @@
 
 
-#include "ternary.h"
+#include "freenergy.h"
 #include "config.h" /*+ For inline +*/
 #include <math.h> /*+ For log, various other functions +*/
@@ -435 +435 @@
   int free_energies It returns zero or an error code.
 
-  ternary_point *points Ternary point structures holding the compositions where
+  energy_point *points Energy point structures holding the compositions where
   we calculate the free energy, and the free energy field where we put the
   return values.
@@ -460 +460 @@
 
 int free_energies
-(ternary_point *points, int n, double T, double P, energy_params *eparams,
+(energy_point *points, int n, double T, double P, energy_params *eparams,
  int efunc, free_energy_flags flags)
 {

trunk/matml/src/ternary/geomview.c

@@ -79 +79 @@
   int npoints Number of points.
 
-  ternary_point *points Point coordinates and colors.
+  visual_point *points Point coordinates and colors.
 
   int ntriangles Number of triangles.
@@ -88 +88 @@
 int GeomviewDisplayTriangleCOFF
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, int ntriangles, int *vertices)
+ int npoints, visual_point *points, int ntriangles, int *vertices)
 {
   int i;
@@ -128 +128 @@
   int npoints Number of points.
 
-  ternary_point *points Point coordinates and colors.
+  visual_point *points Point coordinates and colors.
 
   int numfacets Number of facet objects to display.
@@ -139 +139 @@
 int GeomviewDisplayFacets
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, int numfacets, hull_facet *facets,
+ int npoints, visual_point *points, int numfacets, hull_facet *facets,
  facet_type type)
 {
@@ -189 +189 @@
   int npoints Number of points.
 
-  ternary_point *points Point coordinates and colors.
+  visual_point *points Point coordinates and colors.
 
   phase_boundary theboundary Phase boundary object to display.
@@ -196 +196 @@
 int GeomviewDisplayPhaseBoundary
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, phase_boundary *theboundary)
+ int npoints, visual_point *points, phase_boundary *theboundary)
 {
   int i, c=theboundary->compos;

trunk/matml/src/ternary/Makefile.am

@@ -9 +9 @@
 AUTOMAKE_OPTIONS=1.5
 
+lib_LTLIBRARIES = libfreenergy.la libternary.la
+
+libfreenergy_la_SOURCES = freenergy.c
+libfreenergy_la_CFLAGS  = -std=c99 #-DDEBUG
+libfreenergy_la_LDFLAGS = -version-info 0:0:0
+
+libternary_la_SOURCES = geomview.c qhull.c book.c spinodal.c
+libternary_la_LIBADD  = -lm -lqhull libfreenergy.la
+libternary_la_CFLAGS  = -std=c99 -DGEOMVIEW=\"@GEOMVIEW@\" #-DDEBUG
+libternary_la_LDFLAGS = -version-info 0:0:0
+
+include_HEADERS = freenergy.h ternary.h
+
 bin_PROGRAMS    = ternary square
-noinst_HEADERS  = ternary.h
-ternary_SOURCES = ternary.c geomview.c qhull.c book.c freenergy.c spinodal.c
-ternary_CFLAGS  = -std=c99 -DGEOMVIEW=\"@GEOMVIEW@\" #-DDEBUG
-ternary_LDADD   = -lm -lqhull
-square_SOURCES = square.c geomview.c qhull.c book.c freenergy.c spinodal.c
-square_CFLAGS  = -std=c99 -DGEOMVIEW=\"@GEOMVIEW@\" #-DDEBUG
-square_LDADD   = -lm -lqhull
+ternary_SOURCES = ternary.c
+ternary_CFLAGS  = #-DDEBUG
+ternary_LDADD   = libternary.la libfreenergy.la
+square_SOURCES = square.c
+square_CFLAGS  = #-DDEBUG
+square_LDADD   = libternary.la libfreenergy.la
 
 EXTRA_DIST = autogen.sh macros/autogen.sh macros/cxref-latex.m4 \
@@ -25 +37 @@
 if HAVE_CXREF
 BUILT_TEXFILES = $(top_srcdir)/ternary.c.tex \
+                 $(top_srcdir)/square.c.tex \
                  $(top_srcdir)/ternary.h.tex \
                  $(top_srcdir)/spinodal.c.tex \
@@ -30 +43 @@
                  $(top_srcdir)/qhull.c.tex \
                  $(top_srcdir)/book.c.tex \
+                 $(top_srcdir)/freenergy.h.tex \
                  $(top_srcdir)/freenergy.c.tex \
                  config.h.tex

trunk/matml/src/ternary/qhull.c

@@ -7 +7 @@
 
 
+#include "freenergy.h"
 #include "ternary.h"
 #include <qhull/qset.h>
@@ -102 +103 @@
   int dim Number of dimensions (should always be 3 for ternary).
 
-  ternary_point *points Point information.
+  energy_point *points Point information.
 
   int numpoints Number of points of which to calculate the convex hull.
@@ -121 +122 @@
 
 int hullCalculate
-(int dim, ternary_point *points, int numpoints, energy_params *eparams,
+(int dim, energy_point *points, int numpoints, energy_params *eparams,
  int nparams, double T, double P, hull_facet **facets, int *numfacets)
 {
@@ -312 +313 @@
   int i;
   double distance = 1., G, slope2, slope3;
-  ternary_point current;
+  energy_point current;
 
   if (globaldims != 3)
@@ -350 +351 @@
 
   while (distance >= newton_tolerance * newton_tolerance &&
-         current.C2>=0. && current.C3>=0. && current.C2<=1. && current.C3<=1.)
+         current.C2>=0. && current.C3>=0. &&
+         current.C2<=1. && current.C3<=1.)
     {
       double dC2, dC3, det;
@@ -389 +391 @@
   int hullRefine It returns zero or an error code.
 
-  ternary_point **points Pointer to point list, which will be modified to
+  energy_point **points Pointer to point list, which will be modified to
   accommodate the new points created in the refining process.
 
@@ -422 +424 @@
 
 int hullRefine
-(ternary_point **points, int *numpoints, hull_facet **facets, int *numfacets,
+(energy_point **points, int *numpoints, hull_facet **facets, int *numfacets,
  energy_params *eparams, int nparams, double T,double P,
  double newton_tolerance, double vertex_tolerance, int one_phase_refine)
 {
-  ternary_point *newpoints;
+  energy_point *newpoints;
   int newpointsize=100, numnewpoints=0, i;
 
-  if (!(newpoints = malloc (100*sizeof(ternary_point))))
+  if (!(newpoints = (energy_point *) malloc (100*sizeof(energy_point))))
     return 1;
 
@@ -490 +492 @@
 #endif
           if (numnewpoints >= newpointsize)
-            newpoints = realloc
-              (newpoints, (newpointsize+=100)*sizeof(ternary_point));
+            newpoints = (energy_point *) realloc
+              (newpoints, (newpointsize+=100)*sizeof(energy_point));
           newpoints [numnewpoints].C2 = centroid[0];
           newpoints [numnewpoints].C3 = centroid[1];
@@ -588 +590 @@
 #endif
               if (numnewpoints >= newpointsize)
-                newpoints = realloc
-                  (newpoints, (newpointsize+=100)*sizeof(ternary_point));
+                newpoints = (energy_point *) realloc
+                  (newpoints, (newpointsize+=100)*sizeof(energy_point));
               newpoints [numnewpoints].C2 = newcorners[0];
               newpoints [numnewpoints].C3 = newcorners[1];
@@ -609 +611 @@
 #endif
               if (numnewpoints >= newpointsize)
-                newpoints = realloc
-                  (newpoints, (newpointsize+=100)*sizeof(ternary_point));
+                newpoints = (energy_point *) realloc
+                  (newpoints, (newpointsize+=100)*sizeof(energy_point));
               newpoints [numnewpoints].C2 = newcorners[3];
               newpoints [numnewpoints].C3 = newcorners[4];
@@ -634 +636 @@
 #endif
               if (numnewpoints >= newpointsize)
-                newpoints = realloc
-                  (newpoints, (newpointsize+=100)*sizeof(ternary_point));
+                newpoints = (energy_point *) realloc
+                  (newpoints, (newpointsize+=100)*sizeof(energy_point));
               newpoints [numnewpoints].C2 = newcorners[6];
               newpoints [numnewpoints].C3 = newcorners[7];
@@ -653 +655 @@
 #endif
 
-  if (!(*points = realloc
-        (*points, (*numpoints + numnewpoints) * sizeof (ternary_point))))
+  if (!(*points = (energy_point *) realloc
+        (*points, (*numpoints + numnewpoints) * sizeof (energy_point))))
     return 1;
           
@@ -707 +709 @@
   int hullReturnPhaseBoundaries It returns zero or an error code.
 
-  ternary_point *points Array of points in the ternary system.
+  energy_point *points Array of points in the ternary system.
 
   int n_points Number of points.
@@ -729 +731 @@
 
 int hullReturnPhaseBoundaries
-(ternary_point *points, int n_points, hull_facet *facets, int n_facets,
+(energy_point *points, int n_points, hull_facet *facets, int n_facets,
  energy_params *eparams, double T, double P,
  phase_boundary **boundaries, int *n_bounds)

trunk/matml/src/ternary/spinodal.c

@@ -6 +6 @@
 
 
+#include "freenergy.h"
 #include "ternary.h"
 #include <stdlib.h>  /*+ For malloc() and free() +*/
@@ -42 +43 @@
   int calc_spinodal It returns zero or an error code.
 
-  ternary_point **points Points in the array, this function adds points to the
+  energy_point **points Points in the array, this function adds points to the
   end of the array.
 
@@ -59 +60 @@
 
 int calc_spinodal
-(ternary_point **points, int *n_points, int *triangle_indices, int n_triangles,
+(energy_point **points, int *n_points, int *triangle_indices, int n_triangles,
  double T, double P, energy_params *eparams, int efunc,
  phase_boundary *spinreturn)
@@ -65 +66 @@
   int tri, numsegs=0, thesegs_size=20, *thesegs, numnewpoints=0,
     thenewpoints_size=20;
-  ternary_point *thenewpoints;
+  energy_point *thenewpoints;
 
   //printf ("Calculating free energies at all points\n");
@@ -76 +77 @@
   printf ("Allocating new points and segments arrays\n");
 #endif
-  if (!(thenewpoints = malloc (thenewpoints_size * sizeof (ternary_point))))
+  if (!(thenewpoints = (energy_point *) malloc
+        (thenewpoints_size * sizeof (energy_point))))
     return 1;
   if (!(thesegs = (int *) malloc (2 * thesegs_size * sizeof (int))))
@@ -89 +91 @@
       double spin[3], A2=-1, A3, B2=-1, B3;
       int total;
-      ternary_point P0=(*points)[triangle_indices[3*tri]],
+      energy_point P0=(*points)[triangle_indices[3*tri]],
         P1=(*points)[triangle_indices[3*tri+1]],
         P2=(*points)[triangle_indices[3*tri+2]];
@@ -268 +270 @@
         {
           if (++numnewpoints >= thenewpoints_size)
-            if (!(thenewpoints = (ternary_point *)
-                  malloc ((thenewpoints_size*=2)*sizeof (ternary_point))))
+            if (!(thenewpoints = (energy_point *)
+                  malloc ((thenewpoints_size*=2)*sizeof (energy_point))))
               { free (thenewpoints); free (thesegs); return 1; }
           thenewpoints [numnewpoints-1] .C2 = A2;
@@ -284 +286 @@
         {
           if (++numnewpoints >= thenewpoints_size)
-            if (!(thenewpoints = (ternary_point *)
-                  malloc ((thenewpoints_size*=2)*sizeof (ternary_point))))
+            if (!(thenewpoints = (energy_point *)
+                  malloc ((thenewpoints_size*=2)*sizeof (energy_point))))
               { free (thenewpoints); free (thesegs); return 1; }
           thenewpoints [numnewpoints-1] .C2 = B2;
@@ -312 +314 @@
 
   // Add new points to old points array
-  if (!(*points = realloc
-        (*points, (*n_points+numnewpoints)*sizeof (ternary_point))))
+  if (!(*points = (energy_point *) realloc
+        (*points, (*n_points+numnewpoints)*sizeof (energy_point))))
     { free (thenewpoints); free (thesegs); return 1; }
   for (tri=0; tri<numnewpoints; tri++)
@@ -331 +333 @@
 #endif
       spinreturn->n_edges = numsegs;
-      if (!(spinreturn->phases = malloc (sizeof (int))))
+      if (!(spinreturn->phases = (int *) malloc (sizeof (int))))
         { free (thenewpoints); free (thesegs); return 1; }
       *(spinreturn->phases) = efunc;

trunk/matml/src/ternary/square.c

@@ -33 +33 @@
   FILE *pfd = NULL;
   double T=0.7, P=1.;
-  ternary_point *points;
+  energy_point *points;
+  visual_point *tpoints=NULL;
   phase_boundary solid_spinodal, *allbounds=NULL;
   /* eparams: R,T0, G1@T0,G2,G3, C1,C2,C3, M1,M2,M3, O12,O13,O23,O123 */
@@ -53 +54 @@
 
   /* Allocate array pairs for ternary points and triangle vertices */
-  if (!(points = malloc (numpoints * sizeof (ternary_point))))
+  if (!(points = (energy_point *) malloc (numpoints * sizeof (energy_point))))
     { printf ("Cannot allocate point coordinate array\n"); exit (1); }
-  if (!(verts = malloc (res2*res3 * 12 * sizeof (int))))
+  if (!(verts = (int *) malloc (res2*res3 * 12 * sizeof (int))))
     { printf ("Cannot allocate triangle vertex array\n"); exit (1); }
 
@@ -84 +85 @@
 
   /* Calculate the spinodal: first with no return, then more accurate with one
+     Problem: what's the spinodal criterion in a square array?
   if (i=calc_spinodal (&points, &numpoints, verts, res2*res3*4, T, P,
                        allparams, 0, &solid_spinodal))
@@ -90 +92 @@
 
   /* Scale all free energy values */
-  if (i=scale_energy_array (numpoints, points, 0., -.2, 1., .9, NULL, NULL, 1))
+  if (i=energy_to_visual_array (numpoints, points, &tpoints, 0., -.2, 1., .9,
+                                NULL, NULL, 1))
     { printf ("main: Error %d in scale_triangle_array\n", i); exit (i); }
 
@@ -96 +99 @@
   if (i=GeomviewDisplayTriangleCOFF
       (pfd, "Ternary Free Energy", "tfe", "shading smooth",
-       numpoints, points, res2*res3*4, verts))
+       numpoints, tpoints, res2*res3*4, verts))
     { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
-
-  /*
-  if (i=GeomviewDisplayPhaseBoundary
-      (pfd, "Spinodal", "spin", "-face +edge",
-       numpoints, points, &solid_spinodal))
-    { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
-  */
 
   /* Calculate and refine the lower convex hull */
@@ -122 +118 @@
 
   /* Scale and display everything */
-  if (i=scale_energy_array (numpoints, points, 0., -.2, 1., .9, NULL, NULL, 1))
+  if (i=energy_to_visual_array (numpoints, points, &tpoints, 0., -.2, 1., .9,
+                                NULL, NULL, 1))
     { printf ("main: Error %d in scale_triangle_array\n", i); exit (i); }
-
-  /* No longer needed
-  if (i=GeomviewDisplayFacets
-      (pfd, "Binodal and Tie Lines", "ech", "-face +edge",
-       numpoints, points, hullnumfacets, hullfacets, TWO_PHASE))
-    { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
-  */
 
   if (i=hullReturnPhaseBoundaries (points,numpoints, hullfacets,hullnumfacets,
@@ -151 +141 @@
 
       if (ret=GeomviewDisplayPhaseBoundary
-          (pfd, name, label, "-face +edge", numpoints, points, allbounds+i))
+          (pfd, name, label, "-face +edge", numpoints, tpoints, allbounds+i))
         { printf ("main: Error %d in Geomview Display\n", ret); exit (ret); }
     }

trunk/matml/src/ternary/ternary.c

@@ -33 +33 @@
   FILE *pfd = NULL;
   double T=0.7, P=1.;
-  ternary_point *points;
+  energy_point *points;
+  visual_point *tpoints=NULL;
   phase_boundary solid_spinodal, *allbounds=NULL;
   /* eparams: name, R,T0, G1@T0,G2,G3, C1,C2,C3, S1,S2,S3,S4,S5,
@@ -49 +50 @@
 
   /* Allocate array pairs for ternary points and triangle vertices */
-  if (!(points = malloc (numpoints * sizeof (ternary_point))))
+  if (!(points = (energy_point *) malloc (numpoints * sizeof (energy_point))))
     { printf ("Cannot allocate point coordinate array\n"); exit (1); }
-  if (!(verts = malloc (loop_max*loop_max * 6 * sizeof (int))))
+  if (!(verts = (int *) malloc (loop_max*loop_max * 6 * sizeof (int))))
     { printf ("Cannot allocate triangle vertex array\n"); exit (1); }
 
@@ -85 +86 @@
 
   /* Scale all free energy values */
-  if (i=scale_energy_array (numpoints, points, 0., -.2, 1., .9, NULL, NULL, 0))
+  if (i=energy_to_visual_array (numpoints, points, &tpoints, 0., -.2, 1., .9,
+                                NULL, NULL, 0))
     { printf ("main: Error %d in scale_triangle_array\n", i); exit (i); }
 
@@ -91 +93 @@
   if (i=GeomviewDisplayTriangleCOFF
       (pfd, "Ternary Free Energy", "tfe", "shading smooth",
-       numpoints, points, loop_max*loop_max*2, verts))
+       numpoints, tpoints, loop_max*loop_max*2, verts))
     { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
 
   if (i=GeomviewDisplayPhaseBoundary
       (pfd, "Spinodal", "spin", "-face +edge",
-       numpoints, points, &solid_spinodal))
+       numpoints, tpoints, &solid_spinodal))
     { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
 
@@ -115 +117 @@
 
   /* Scale and display everything */
-  if (i=scale_energy_array (numpoints, points, 0., -.2, 1., .9, NULL, NULL, 0))
+  if (i=energy_to_visual_array (numpoints, points, &tpoints, 0., -.2, 1., .9,
+                                NULL, NULL, 0))
     { printf ("main: Error %d in scale_triangle_array\n", i); exit (i); }
-
-  /* No longer needed
-  if (i=GeomviewDisplayFacets
-      (pfd, "Binodal and Tie Lines", "ech", "-face +edge",
-       numpoints, points, hullnumfacets, hullfacets, TWO_PHASE))
-    { printf ("main: Error %d in Geomview Display\n", i); exit (i); }
-  */
 
   if (i=hullReturnPhaseBoundaries (points,numpoints, hullfacets,hullnumfacets,
@@ -144 +140 @@
 
       if (ret=GeomviewDisplayPhaseBoundary
-          (pfd, name, label, "-face +edge", numpoints, points, allbounds+i))
+          (pfd, name, label, "-face +edge", numpoints, tpoints, allbounds+i))
         { printf ("main: Error %d in Geomview Display\n", ret); exit (ret); }
     }

trunk/matml/src/ternary/ternary.h

@@ -2 +2 @@
   $Header$
 
-  This is the header file.
+  This is the ternary.h header file with structures and prototypes for
+  calculating ternary phase diagrams.
   ***************************************/
 
@@ -10 +11 @@
 
 #include <stdio.h>
+#include "freenergy.h"
 
 typedef struct {
-  double C2;     /*+ Ternary concentration variable
-                   +latex+$C_2$ ($0\rightarrow1$)
-                   -latex-2 (0 -> 1)
-                   +*/
-  double C3;     /*+ Ternary concentration variable
-                   +latex+$C_3$ ($0\rightarrow1-C_2$)
-                   -latex-3 (0 -> 1-C2)
-                   +*/
-  int efunc;     /*+ Index indicating which energy surface (phase) this point
-                   is on; a negative value indicates that energy derivatives
-                   are not yet calculated. +*/
-  double T;      /*+ Temperature at which the free energies and derivatives
-                   are valid +*/
-  double P;      /*+ Pressure at which the free energies and derivatives are
-                   valid, -1 if not valid +*/
-  double G;      /*+ Free energy +*/
-  double G2;     /*+ Free energy derivative
-                   +latex+$\partial G/\partial C_2$
-                   -latex-dG/dC2
-                   +*/
-  double G3;     /*+ Free energy derivative
-                   +latex+$\partial G/\partial C_3$
-                   -latex-dG/dC3
-                   +*/
-  double G22;    /*+ Free energy second derivative
-                   +latex+$\partial^2 G/\partial C_2^2$
-                   +html+ d<sup>2</sup>G/dC2<sup>2</sup>
-                   +*/
-  double G33;    /*+ Free energy second derivative
-                   +latex+$\partial^2 G/\partial C_3^2$
-                   +html+ d<sup>2</sup>G/dC3<sup>2</sup>
-                   +*/
-  double G23;    /*+ Free energy mixed partial second derivative
-                   +latex+$\partial^2 G/\partial C_2\partial C_3$
-                   +html+ d<sup>2</sup>G/dC2dC3
-                   +*/
   double x;      /*+ Visualization
                    +latex+$x$-coordinate ($0\rightarrow1$)
@@ -75 +41 @@
                    -latex-(0 -> 1)
                    +*/
-} ternary_point; /*+ Ternary "point" structure. +*/
-
-typedef struct {
-  char *name;
-  double C2; /*+ Composition at center of this Gaussian term +*/
-  double C3; /*+ Composition at center of this Gaussian term +*/
-  double G;  /*+ Free energy addition at the Gaussian center +*/
-  double w;  /*+ Width of the Gaussian distribution +*/
-  double n;  /*+ Exponent of the distribution, currently fixed at 2 +*/
-} energy_gaussian;
-
-typedef struct {
-  char *name;      /*+ Name of this phase +*/
-  double R;        /*+ Ideal gas law constant +*/
-  double T0;       /*+ Reference temperature +*/
-  double G1_T0;    /*+ Free energy of species 1 at reference temperature +*/
-  double G2_T0;    /*+ Free energy of species 2 at reference temperature +*/
-  double G3_T0;    /*+ Free energy of species 3 at reference temperature +*/
-  double G1_C;     /*+ Species 1 heat capacity +*/
-  double G2_C;     /*+ Species 2 heat capacity +*/
-  double G3_C;     /*+ Species 3 heat capacity +*/
-  double S1;       /*+ Species 1 entropy coefficient (inverse Flory-Huggins relative molar volume) +*/
-  double S2;       /*+ Species 2 entropy coefficient +*/
-  double S3;       /*+ Species 3 entropy coefficient +*/
-  double S4;       /*+ 1-C2 entropy coefficient for square system +*/
-  double S5;       /*+ 1-C3 entropy coefficient for square system +*/
-  double Omega12;  /*+ Species 1-2 regular solution interaction parameter +*/
-  double Omega13;  /*+ Species 1-3 regular solution interaction parameter +*/
-  double Omega23;  /*+ Species 2-3 regular solution interaction parameter +*/
-  double Omega123; /*+ Species 1-2-3 regular solution interaction parameter +*/
-  double Omega234; /*+ Species 2-3-4 regular solution interaction parameter +*/
-  double Omega235; /*+ Species 2-3-5 regular solution interaction parameter +*/
-  double Omega245; /*+ Species 2-4-5 regular solution interaction parameter +*/
-  double Omega345; /*+ Species 3-4-5 regular solution interaction parameter +*/
-  double Omega2345;/*+ Species 2-3-4-5 regular solution interaction parameter+*/
-  energy_gaussian *gauss; /*+ Array of Gaussian energy parameter structs +*/
-  int n_gauss;     /*+ Number of Gaussian energy terms +*/
-} energy_params;   /*+ Free energy parameters structure +*/
-
-typedef enum {
-  STANDARD=0,
-  WITH_DERIVATIVES=1,
-  NO_INFINITY=2,
-  WITH_DERIVATIVES_NO_INFINITY=3
-} free_energy_flags;
+} visual_point; /*+ Ternary "point" structure. +*/
 
 typedef enum {
@@ -130 +52 @@
 typedef struct {
   int compos;   /*+ Number of points in a boundary edge +*/
-  int *edges;   /*+ Edges with indices of ternary_points comprising a phase
+  int *edges;   /*+ Edges with indices of energy_points comprising a phase
                   boundary; each edge has compos points +*/
   int n_edges;  /*+ Number of edges comprising this boundary +*/
@@ -140 +62 @@
                     boundary between a single phase and a multi-phase region,
                     or the tie lines/simplices across a multi-phase region +*/
-
-#define TERNARY_INFINITY HUGE_VALF
-#define TERNARY_NEGATIVE_INFINITY (-HUGE_VALF)
 
 typedef enum {
@@ -157 +76 @@
 } hull_facet;
 
-/* From freenergy.c */
-double free_energy
-(double C2, double C3, double T, double P, energy_params *eparams);
-int free_energies
-(ternary_point *points, int n, double T, double P, energy_params *eparams,
- int efunc, free_energy_flags flags);
-
 /* From spinodal.c */
 int calc_spinodal
-(ternary_point **points, int *n_points, int *triangle_indices, int n_triangles,
+(energy_point **points, int *n_points, int *triangle_indices, int n_triangles,
  double T, double P, energy_params *eparams, int efunc,
  phase_boundary *spinreturn);
@@ -174 +86 @@
 int GeomviewDisplayTriangleCOFF
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, int ntriangles, int *vertices);
+ int npoints, visual_point *points, int ntriangles, int *vertices);
 int GeomviewDisplayFacets
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, int numfacets, hull_facet *facets,
+ int npoints, visual_point *points, int numfacets, hull_facet *facets,
  facet_type type);
 int GeomviewDisplayPhaseBoundary
 (FILE *geompipe, char *name, char *label, char *appearance,
- int npoints, ternary_point *points, phase_boundary *theboundary);
+ int npoints, visual_point *points, phase_boundary *theboundary);
 int GeomviewEnd (FILE **geompipe);
 
 /* From qhull.c */
 int hullCalculate
-(int dim, ternary_point *points, int numpoints, energy_params *eparams,
+(int dim, energy_point *points, int numpoints, energy_params *eparams,
  int nparams, double T, double P, hull_facet **facets, int *numfacets);
 int hullRefine
-(ternary_point **points, int *numpoints, hull_facet **facets, int *numfacets,
+(energy_point **points, int *numpoints, hull_facet **facets, int *numfacets,
  energy_params *eparams, int nparams, double T,double P,
  double newton_tolerance, double vertex_tolerance, int one_phase_refine);
 int hullReturnPhaseBoundaries
-(ternary_point *points, int n_points, hull_facet *facets, int n_facets,
+(energy_point *points, int n_points, hull_facet *facets, int n_facets,
  energy_params *eparams, double T, double P,
  phase_boundary **boundaries, int *n_bounds);
@@ -200 +112 @@
 /* From book.c */
 int init_triangle_array
-(int resolution, ternary_point *points, int efunc,
+(int resolution, energy_point *points, int efunc,
  double A2, double A3, double B2, double B3, double C2, double C3);
 int init_rectangle_array
-(int res2, int res3, ternary_point *points, int efunc,
+(int res2, int res3, energy_point *points, int efunc,
  double A2, double A3, double B2, double B3);
 int init_triangle_vertices (int resolution, int *vertex_array, int offset);
 int init_rectangle_vertices (int res2, int res3, int *vertex_array, int offset);
-int scale_energy_array (int num_points, ternary_point *points,
-                        double y0, double G0, double y1, double G1,
-                        double *Gmin, double *Gmax, int square);
+int energy_to_visual_array
+(int num_points, energy_point *epoints, visual_point **tpoints, double y0,
+ double G0, double y1, double G1, double *Gmin, double *Gmax, int square);
 
 #endif /* TERNARY_H */

trunk/matml/src/ternary/ternary.tex.in

@@ -46 +46 @@
 \section{How it works}
 
-The file {\tt ternary.c} (documentation in appendix \ref{file_ternary.c})
-contains {\tt main()}, which is a simple demonstration of the ternary
-calculation and visualization functions.  It does the following:
+The files {\tt ternary.c} and {\tt square.c} (documentation in appendices
+\ref{file_ternary.c} and \ref{file_square.c}) contain sample front-ends to the
+free energy and phase diagram libraries.  They do the following:
 \begin{itemize}
-\item Calls functions in {\tt book.c} (appendix \ref{file_book.c}) to create a
-  triangular array of points and set of triangles connecting them.
-\item Calculates the free energy function on the triangle verticies using
+\item Call functions in {\tt book.c} (appendix \ref{file_book.c}) to create a
+  triangular or square array of points and set of triangles connecting them.
+\item Calculate the free energy function on the triangle vertices using
   {\tt free\_energies()} in {\tt freenergy.c} (appendix
   \ref{file_freenergy.c}).
-\item Displays the free energy using the functions in {\tt geomview.c}
+\item Display the free energy using the functions in {\tt geomview.c}
   (appendix \ref{file_geomview.c}), which forks and controls a Geomview
   process.
-\item Calculates the spinodal of this free energy function using {\tt
+\item Calculate the spinodal of this free energy function using {\tt
     calc\_spinodal} in {\tt spinodal.c} (appendix \ref{file_spinodal.c}), and
   displays it on the free energy function.
-\item Uses Qhull calls in {\tt qhull.c} (appendix \ref{file_qhull.c}) to
+\item Use Qhull calls in {\tt qhull.c} (appendix \ref{file_qhull.c}) to
   calculate the convex hull of the free energy function.
-\item Refines the multi-phase triangles using Newton iteration ({\tt
+\item Refine the multi-phase triangles using Newton iteration ({\tt
     hullRefine} in {\tt qhull.c}) to find the lowest-energy point starting at
   each corner, smoothing the phase boundaries considerably.
-\item Adds the two-phase region facet outlines to the Geomview display.
+\item Add the two-phase region facet outlines to the Geomview display.
 \end{itemize}
 
@@ -105 +105 @@
   spinodal.c} calculates the shape of the spinodal curve.
 
+All functions can now accommodate ``rectangle'' ternary spaces, though it's not
+clear what a spinodal means in this context.  The {\tt square} front-end
+demonstrates this functionality.
+
+Finally, a major reorganization splits the bulk of the code into two libraries,
+with front-end codes {\tt ternary} and {\tt square} calling into them.  One
+library calculates free energies and derivatives, and can be re-used in phase
+field codes.  The second calculates phase diagrams.  The APIs are relatively
+immature and bound to change.
+
 \bibliographystyle{unsrturl}
 \bibliography{ternary}