Commits

John Wise committed efd0b80 Merge

Comments (0)

Files changed (7)

src/enzo/FluxFix_Grid_CorrectForRefinedFluxes.C

  
 int FindField(int f, int farray[], int n);
 int CosmologyComputeExpansionFactor(FLOAT time, FLOAT *a, FLOAT *dadt);
+int MakeFieldConservative(field_type field); 
  
 int grid::CorrectForRefinedFluxes(fluxes *InitialFluxes,
 				  fluxes *RefinedFluxes,
 	  if (HydroMethod != Zeus_Hydro) {
 
 	    for (field = 0; field < NumberOfBaryonFields; field++) 
-	      if (FieldTypeNoInterpolate(FieldType[field]) == FALSE &&
-		  FieldTypeIsDensity(FieldType[field]) == FALSE &&
-		  FieldTypeIsRadiation(FieldType[field]) == FALSE &&
-		  FieldType[field] != Bfield1 &&
-		  FieldType[field] != Bfield2 && FieldType[field] != Bfield3 &&
-		  FieldType[field] != PhiField && 
-		  FieldType[field] != DrivingField1 &&
-		  FieldType[field] != DrivingField2 &&
-		  FieldType[field] != DrivingField3 &&
-		  FieldType[field] != GravPotential &&
-		  FieldType[field] != DebugField &&
-		(RadiativeCooling == 0 || (FieldType[field] != TotalEnergy && 
-					   FieldType[field] != InternalEnergy))) {
+	      if ( MakeFieldConservative(FieldType[field]) ) {
 		for (k = Start[2]; k <= End[2]; k++) {
 		  for (j = Start[1]; j <= End[1]; j++) {
 		    index = (k*GridDimension[1] + j)*GridDimension[0] + Start[0];
 	
 	  if (HydroMethod != Zeus_Hydro)
 	    for (field = 0; field < NumberOfBaryonFields; field++)
-	      if (FieldTypeIsDensity(FieldType[field]) == FALSE &&
-		  FieldTypeNoInterpolate(FieldType[field]) == FALSE &&
-		  FieldTypeIsRadiation(FieldType[field]) == FALSE &&
-		  (RadiativeCooling == 0 || (FieldType[field] != TotalEnergy &&
-					     FieldType[field] != InternalEnergy)) && 
-		  FieldType[field] != Bfield1 &&
-		  FieldType[field] != Bfield2 && FieldType[field] != Bfield3 &&
-		  FieldType[field] != PhiField &&
-		  FieldType[field] != DrivingField1 &&
-		  FieldType[field] != DrivingField2 &&
-		  FieldType[field] != DrivingField3 &&
-		  FieldType[field] != GravPotential)
+	      if ( MakeFieldConservative(FieldType[field]))
 		for (k = Start[2]; k <= End[2]; k++)
 		  for (j = Start[1]; j <= End[1]; j++) {
 		    index = (k*GridDimension[1] + j)*GridDimension[0] + Start[0];

src/enzo/Grid_CorrectForRefinedFluxes.C

  
 	if (HydroMethod != Zeus_Hydro)
 	  for (field = 0; field < NumberOfBaryonFields; field++)
-	    if (FieldTypeIsDensity(FieldType[field]) == FALSE &&
-		FieldTypeIsRadiation(FieldType[field]) == FALSE &&
-		FieldType[field] != Bfield1 &&
-		FieldType[field] != Bfield2 && FieldType[field] != Bfield3 &&
-		FieldType[field] != PhiField &&
-		FieldType[field] != DrivingField1 &&
-		FieldType[field] != DrivingField2 &&
-		FieldType[field] != DrivingField3 &&
-		FieldType[field] != GravPotential)  
+	    if ( MakeFieldConservative(FieldType[field]) )
+                
 	      //		(RadiativeCooling == 0 || (FieldType[field] != TotalEnergy &&
 	      //	 			 FieldType[field] != InternalEnergy)))
 	      for (k = Start[2]; k <= End[2]; k++)
  
 	if (HydroMethod != Zeus_Hydro)
 	  for (field = 0; field < NumberOfBaryonFields; field++)
-	    if (FieldTypeIsDensity(FieldType[field]) == FALSE &&
-		FieldTypeIsRadiation(FieldType[field]) == FALSE &&
-		FieldType[field] != Bfield1 &&
-		FieldType[field] != Bfield2 && FieldType[field] != Bfield3 &&
-		FieldType[field] != PhiField &&
-		FieldType[field] != DrivingField1 &&
-		FieldType[field] != DrivingField2 &&
-		FieldType[field] != DrivingField3 &&
-		FieldType[field] != GravPotential)
+	    if ( MakeFieldConservative(FieldType[field]))
+        
 	      for (k = Start[2]; k <= End[2]; k++)
 		for (j = Start[1]; j <= End[1]; j++) {
 		  index = (k*GridDimension[1] + j)*GridDimension[0] + Start[0];

src/enzo/Grid_InterpolateBoundaryFromParent.C

 	       int *ivel_flag, int *irefine);
  
 /* InterpolateBoundaryFromParent function */
- 
+int MakeFieldConservative(field_type field); 
 int grid::InterpolateBoundaryFromParent(grid *ParentGrid)
 {
  
  
     if (ConservativeInterpolation)
       for (field = 0; field < NumberOfBaryonFields; field++)
-	if (FieldTypeIsDensity(FieldType[field]) == FALSE &&
-	    FieldType[field] != Bfield1 &&
-	    FieldType[field] != Bfield2 &&
-	    FieldType[field] != Bfield3 &&
-	    FieldType[field] != PhiField &&
-	    FieldType[field] != DrivingField1 &&
-	    FieldType[field] != DrivingField2 &&
-	    FieldType[field] != DrivingField3 &&
-	    FieldType[field] != DebugField &&
-	    FieldType[field] != GravPotential
-	    ) {
+	if (MakeFieldConservative(FieldType[field])) {
 	  FORTRAN_NAME(mult3d)(ParentTemp[densfield], ParentTemp[field],
 			       &ParentTempSize, &One, &One,
 			       &ParentTempSize, &One, &One,
          variables (skipping density). */
  
       if (ConservativeInterpolation)
-	if (FieldTypeIsDensity(FieldType[field]) == FALSE  &&
-	    FieldType[field] != Bfield1 &&
-	    FieldType[field] != Bfield2 &&
-	    FieldType[field] != Bfield3 &&
-	    FieldType[field] != PhiField &&
-	    FieldType[field] != DrivingField1 &&
-	    FieldType[field] != DrivingField2 &&
-	    FieldType[field] != DrivingField3 &&
-	    FieldType[field] != DebugField &&
-	    FieldType[field] != GravPotential
-	    )
+	if (MakeFieldConservative(FieldType[field]))
 	  FORTRAN_NAME(div3d)(TemporaryDensityField, TemporaryField,
 			      &TempSize, &One, &One,
 			      &TempSize, &One, &One,

src/enzo/Grid_InterpolateFieldValues.C

 	       int *dstart1, int *dstart2, int *dstart3,
 	       int *ivel_flag, int *irefine);
  
+int MakeFieldConservative(field_type field); 
+
 /* InterpolateBoundaryFromParent function */
  
 int grid::InterpolateFieldValues(grid *ParentGrid)
        quantities. */
  
     if (ConservativeInterpolation)
-      for (field = 0; field < NumberOfBaryonFields; field++)
-	if (FieldTypeIsDensity(FieldType[field]) == FALSE &&
-	    FieldType[field] != Bfield1 &&
-	    FieldType[field] != Bfield2 &&
-	    FieldType[field] != Bfield3 &&
-	    FieldType[field] != PhiField &&
-	    FieldType[field] != DrivingField1 &&
-	    FieldType[field] != DrivingField2 &&
-	    FieldType[field] != DrivingField3 
-	    && FieldType[field] != DebugField 
-	    && FieldType[field] != GravPotential
-	    )
+      for (field = 0; field < NumberOfBaryonFields; field++){
+	if (MakeFieldConservative( FieldType[field] ) ){
 	  FORTRAN_NAME(mult3d)(ParentTemp[densfield], ParentTemp[field],
                                &ParentTempSize, &One, &One,
 			       &ParentTempSize, &One, &One,
                                &Zero, &Zero, &Zero, &Zero, &Zero, &Zero);
+    }
+      }
     
     /* Do the interpolation for the density field. */
  
          variables (skipping density). */
  
       if (ConservativeInterpolation)
-	if (FieldTypeIsDensity(FieldType[field]) == FALSE  &&
-	    FieldType[field] != Bfield1 &&
-	    FieldType[field] != Bfield2 &&
-	    FieldType[field] != Bfield3 &&
-	    FieldType[field] != PhiField &&
-	    FieldType[field] != DrivingField1 &&
-	    FieldType[field] != DrivingField2 &&
-	    FieldType[field] != DrivingField3 
-	    && FieldType[field] != DebugField 
-	    &&  FieldType[field] != GravPotential
-	    )
+	if (MakeFieldConservative( FieldType[field] ) ){
 	  FORTRAN_NAME(div3d)(TemporaryDensityField, TemporaryField,
 			      &TempSize, &One, &One,
 			      &TempSize, &One, &One,
 			      &Zero, &Zero, &Zero, &Zero, &Zero, &Zero,
 			      &Zero, &Zero, &Zero, &TempSize, &Zero, &Zero);
+    }
       
       /* Set FieldPointer to either the correct field (density or the one we
 	 just interpolated to). */

src/enzo/Grid_ProjectSolutionToParentGrid.C

 				    int *dstart1, int *dstart2, int *dstart3,
                                     int *rstart1, int *rstart2, int *rstart3,
                                     int *rend1, int *rend2, int *rend3);
- 
- 
+int MakeFieldConservative(field_type field);
 int grid::ProjectSolutionToParentGrid(grid &ParentGrid)
 {
   /* Return if this doesn't involve us. */
  
   if (ProcessorNumber == MyProcessorNumber)
     for (field = 0; field < NumberOfBaryonFields; field++)
-      if (FieldTypeIsDensity(FieldType[field]) == FALSE && 
-	  (FieldTypeNoInterpolate(FieldType[field]) == FALSE) && (
-	  (FieldType[field] < Velocity1 || FieldType[field] > Velocity3)
-          || HydroMethod != Zeus_Hydro) &&
-	  FieldType[field] != Bfield1 &&
-	  FieldType[field] != Bfield2 &&
-	  FieldType[field] != Bfield3 &&
-	  FieldType[field] != PhiField &&
-	  FieldType[field] != DrivingField1 &&
-	  FieldType[field] != DrivingField2 &&
-	  FieldType[field] != DrivingField3 &&
-	  FieldType[field] != GravPotential)
+      if (MakeFieldConservative(FieldType[field]))
       FORTRAN_NAME(mult3d)(BaryonField[DensNum], BaryonField[field],
 			   &Size, &One, &One, &Size, &One, &One,
 			   &Zero, &Zero, &Zero, &Zero, &Zero, &Zero);
   /* Divide all fields by mass to return to original quantity. */
  
   for (field = 0; field < NumberOfBaryonFields; field++)
-    if ((FieldTypeIsDensity(FieldType[field]) == FALSE &&
-	 (FieldTypeNoInterpolate(FieldType[field]) == FALSE) && (
-	  (FieldType[field] < Velocity1 || FieldType[field] > Velocity3)
-          || HydroMethod != Zeus_Hydro) &&
-	FieldType[field] != Bfield1 &&
-	FieldType[field] != Bfield2 &&
-	FieldType[field] != Bfield3 &&
-	FieldType[field] != PhiField &&
-	FieldType[field] != DrivingField1 &&
-	FieldType[field] != DrivingField2 &&
-	FieldType[field] != DrivingField3 &&
-	 FieldType[field] != GravPotential )) {
+    if ( MakeFieldConservative(FieldType[field]) ) {
       if (ProcessorNumber == MyProcessorNumber)
 	FORTRAN_NAME(div3d)(BaryonField[DensNum], BaryonField[field],
 			    &Size, &One, &One, &Size, &One, &One,

src/enzo/Make.config.objects

         MagneticFieldResetter.o \
         mbh_maker.o \
         mcooling.o \
+        MakeFieldConservative.o\
         MemoryAllocationRoutines.o \
 	MemoryPoolRoutines.o \
 	MersenneTwister.o \

src/enzo/MakeFieldConservative.C

+/* 
+ * MakeFieldConservative
+ *
+ * Returns True if field should be multiplied by Density for various operations
+ * that typically expect conservative fields.  This includes interpolation and flux correction.
+ *
+ * David Collins. 2011-05-25
+ */
+#include <stdio.h>
+#include "ErrorExceptions.h"
+#include "macros_and_parameters.h"
+#include "typedefs.h"
+#include "global_data.h"
+#include "Fluxes.h"
+#include "GridList.h"
+#include "ExternalBoundary.h"
+#include "Grid.h"
+
+
+int MakeFieldConservative(field_type field){
+    int MultiplyField = TRUE, counter = 0;
+    //types_to_skip lists types not to be multiplied.
+    //Please leave FieldUndefined as the last element, it terminates the loop.
+    //
+    field_type  types_to_skip[] = {Bfield1, Bfield2, Bfield3, PhiField,
+        DrivingField1, DrivingField2, DrivingField3, GravPotential, DebugField,
+    FieldUndefined};
+
+    if( FieldTypeIsDensity(field) ) MultiplyField = FALSE;
+    if( FieldTypeIsRadiation(field) ) MultiplyField = FALSE;
+    if( FieldTypeNoInterpolate( field ) ) MultiplyField = FALSE;
+    if( (field >= Velocity1 && field <= Velocity3) && HydroMethod == Zeus_Hydro ) {
+        MultiplyField = FALSE;
+    }
+    while( types_to_skip[counter] != FieldUndefined ){
+        if( field == types_to_skip[counter++] ){
+            MultiplyField = FALSE;
+            break;
+        }
+    }
+
+    return MultiplyField;
+}