Merge branch 'devel' into elmerice

.github/workflows/build-macos-homebrew.yaml

@@ -35,10 +35,18 @@ jobs:
         uses: actions/checkout@v4
 
       - name: install dependencies
+
+        # Homebrew's Python conflicts with the Python that comes pre-installed
+        # on the GitHub runners. Some of the dependencies depend on different
+        # versions of Homebrew's Python. Enforce using the ones from Homebrew
+        # to avoid errors on updates.
+        # See: https://github.com/orgs/Homebrew/discussions/3928
+
         # It looks like "gfortran" isn't working correctly unless "gcc" is
         # re-installed.
         run: |
           brew update
+          brew install --overwrite [email protected] [email protected]
           brew reinstall gcc
           brew install \
             cmake libomp openblas open-mpi suitesparse \

.github/workflows/build.yaml

@@ -90,7 +90,8 @@ jobs:
             libqwt-qt5-dev qtscript5-dev libqt5svg5-dev \
             libvtk9-qt-dev libglvnd-dev \
             occt-misc libocct-data-exchange-dev libocct-draw-dev \
-            $([ "${{ matrix.mpi }}" == "with" ] && echo "trilinos-all-dev libptscotch-dev")
+            $([ "${{ matrix.mpi }}" == "with" ] && echo "trilinos-all-dev libptscotch-dev") \
+            $([ "${{ matrix.mpi }}" == "with" ] && echo "libnetcdff-dev")
 
       - name: configure
         # CHOLMOD and rocALUTION require a working OpenMP package. So, disable them for clang.

ElmerGUI/Application/plugins/egconvert.cpp

@@ -5280,17 +5280,24 @@ int LoadGmshInput(struct FemType *data,struct BoundaryType *bound,
   }
 
   if(strstr(line,"$")) {
-    int verno,minorno;
+    int verno,minorno,gmshformat;
     char *cp;
-    
+
     Getrow(line,in,FALSE);
     cp = line;    
     verno = next_int(&cp);
     cp++;
     minorno = next_int(&cp);
 
     if(info) printf("Gmsh version is %d.%d\n",verno,minorno);
-
+
+    cp++;
+    gmshformat = next_int(&cp);
+    if(gmshformat == 1){
+      printf("Error: Gmsh input file is in binary format! Exiting.\n");
+      bigerror("Gmsh input file is in binary format!");
+    }
+
     fclose(in);
 
     if( verno == 4 ) {

cmake/Modules/FindNETCDF.cmake

@@ -64,7 +64,7 @@ IF (NETCDF_INCLUDE_DIR AND NETCDF_LIBRARY AND NETCDFF_LIBRARY)
   UNSET(NETCDF_FAILMSG)
   SET(NETCDFLIB_FOUND TRUE)
   SET(NETCDF_INCLUDE_DIR ${NETCDF_INCLUDE_DIR})
-  SET(NETCDF_LIBRARIES "${NETCDF_LIBRARY};${NETCDFF_LIBRARY}")
+  SET(NETCDF_LIBRARIES "${NETCDFF_LIBRARY};${NETCDF_LIBRARY}")
 ELSE()
   SET(NETCDF_FAILMSG "NETCDF libraries not found.")
 ENDIF()

elmergrid/src/egconvert.c

@@ -5280,17 +5280,24 @@ int LoadGmshInput(struct FemType *data,struct BoundaryType *bound,
   }
 
   if(strstr(line,"$")) {
-    int verno,minorno;
+    int verno,minorno,gmshformat;
     char *cp;
-    
+
     Getrow(line,in,FALSE);
     cp = line;    
     verno = next_int(&cp);
     cp++;
     minorno = next_int(&cp);
 
     if(info) printf("Gmsh version is %d.%d\n",verno,minorno);
-
+
+    cp++;
+    gmshformat = next_int(&cp);
+    if(gmshformat == 1){
+      printf("Error: Gmsh input file is in binary format! Exiting.\n");
+      bigerror("Gmsh input file is in binary format!");
+    }
+
     fclose(in);
 
     if( verno == 4 ) {

elmergrid/src/egextra.c

@@ -820,7 +820,111 @@ int SaveSizeInfo(struct FemType *data,struct BoundaryType *bound,
   return(0);
 }
 
+int MeshPieces(struct FemType *data,int nomesh,int nomeshes,int info)
+{
+  int i,j,k,n;
+  int MinIndex,MaxIndex,NoPieces,NoCorners,Loop,Ready;
+  int *MeshPiece=NULL,*PiecePerm=NULL;
+  // Indexes only needs to hold the max number of dofs per element
+  int Indexes[100];
+
+  if(nomeshes > 1) {
+    printf("Calculate Mesh Pieces in mesh[%d] of [%d] meshes:\n",nomesh,nomeshes);
+  } else {
+    printf("Calculate Mesh Pieces in mesh:\n");
+  }
+
+  n = data->noknots;
+  MeshPiece = Ivector(1,n);
+  for(i=1;i<=n;i++) MeshPiece[i] = 0;
+
+  /*  Only set the piece for the nodes that are used by some element
+      For others the marker will remain zero. */
+  for(i=1; i<=data->noelements; i++) {
+    NoCorners = data->elementtypes[i]%100;
+    for(j=0; j<NoCorners; j++) {
+      MeshPiece[data->topology[i][j]] = 1;
+    }
+  }
+
+  j = 0;
+  for(i=1; i<=n; i++) {
+    if(MeshPiece[i] > 0) {
+      j++;
+      MeshPiece[i] = j;
+    }
+  }
+  if( n > j) {
+    printf("Number of non-body (hanging) nodes in mesh is %d\n", n-j );
+    printf("Consider running ElmerGrid with -autoclean command\n");
+  }
 
+  /* We go through the elements and set all the piece indexes to minimum index
+     until the mesh is unchanged. Thereafter the whole piece will have the minimum index
+     of the piece.  */
+  Ready = FALSE;
+  Loop = 0;
+
+  while(!Ready) {
+    Ready = TRUE;
+    for(i=1; i<=data->noelements; i++) {
+      MaxIndex = 0;
+      MinIndex = n;
+
+      NoCorners = data->elementtypes[i]%100;
+      for(j=0; j<NoCorners; j++) {
+        Indexes[j] = data->topology[i][j];
+      }
+      for(j=0; j<NoCorners; j++) {
+        for(k=0; k<NoCorners; k++) {
+          if(MaxIndex <  MeshPiece[Indexes[k]] )
+            MaxIndex = MeshPiece[Indexes[k]];
+          if(MinIndex >= MeshPiece[Indexes[k]] )
+            MinIndex = MeshPiece[Indexes[k]];
+        }
+        if(MaxIndex > MinIndex) {
+          MeshPiece[Indexes[j]] = MinIndex;
+          Ready = FALSE;
+        }
+      }
+    }
+    Loop++;
+  }
+/*  printf("Mesh coloring loops: %d\n",Loop); */
+
+  MaxIndex = 0;
+  for(i=1; i<=n; i++)
+    if(MeshPiece[i] > MaxIndex) MaxIndex = MeshPiece[i];
+
+  /*  Compute the true number of different pieces  */
+  if(MaxIndex == 1) {
+    NoPieces = 1;
+  } else {
+    NoPieces = 0;
+    PiecePerm = Ivector(1,MaxIndex);
+    for(i=1;i<=MaxIndex;i++) 
+      PiecePerm[i] = 0;    
+
+    for(i=1; i<=n; i++) {
+      j = MeshPiece[i];
+      if( j == 0) continue;
+      if(PiecePerm[j] == 0) {
+        NoPieces++;
+        PiecePerm[j] = NoPieces;
+      }
+    }
+    free_Ivector(PiecePerm,1,MaxIndex);
+  }
+  if(NoPieces == 1) {
+    printf("There is a single piece in the mesh, so the mesh is conforming.\n");
+  } else {
+    printf("Number of separate pieces in mesh is %d\n", NoPieces);
+    printf("The mesh is non-conforming. If not expecting a non-conforming\n");
+    printf("mesh, then refer to the Elmer User Forum for help.\n");
+  }
+  free_Ivector(MeshPiece,1,n);
+  return(0);
+}
 
 int SaveElmerInputFemBem(struct FemType *data,struct BoundaryType *bound,
 			 char *prefix,int decimals,int info)

elmergrid/src/egextra.h

@@ -26,6 +26,7 @@ int SaveSizeInfo(struct FemType *data,struct BoundaryType *bound,
 		 char *prefix,int info);
 int SaveElmerInputFemBem(struct FemType *data,struct BoundaryType *bound,
 			 char *prefix,int decimals, int info);
+int MeshPieces(struct FemType *data,int nomesh,int nomeshes,int info);
 
 void InspectVector(Real *vector,int first,int last,Real *min,
 		   Real *max,int *mini,int *maxi);

elmergrid/src/fempre.c

@@ -784,6 +784,9 @@ int main(int argc, char *argv[])
   if(info) for(k=0;k<nomeshes;k++)
     BoundingBox(&data[k],k,nomeshes,info);
 
+  if(info) for(k=0;k<nomeshes;k++)
+    MeshPieces(&data[k],k,nomeshes,info);
+
   if(eg.nosave) {
     Goodbye();
     return(0);

elmerice/Solvers/CMakeLists.txt

@@ -4,27 +4,27 @@ SET(WITH_ScatteredDataInterpolator FALSE CACHE BOOL "Include ElmerIce ScatteredD
 MARK_AS_ADVANCED(WITH_ScatteredDataInterpolator)
 
 # ---------------------- #
-# -- NetCDF libraries -- #
+# -- netCDF libraries -- #
 #----------------------- #
 MESSAGE(STATUS "------------------------------------------------")
-MESSAGE(STATUS "Elmer/Ice package: Looking for [NetCDF] & [NetCDF Fortran] libraries")
+MESSAGE(STATUS "Elmer/Ice package: Looking for [netCDF] & [netCDF Fortran] libraries")
 
-FIND_PACKAGE(NetCDF)
+FIND_PACKAGE(NETCDF MODULE)
 
-IF(NetCDF_FOUND)
+IF(NETCDF_FOUND)
 
   SET(HAVE_NETCDF TRUE)
   MARK_AS_ADVANCED(HAVE_NETCDF)
-  INCLUDE_DIRECTORIES(${NetCDF_INCLUDE_DIR})
+  INCLUDE_DIRECTORIES(${NETCDF_INCLUDE_DIR})
   ADD_DEFINITIONS(-DHAVE_NETCDF)
 
-  MESSAGE(STATUS "  NetCDF:           " "${NetCDF_FOUND}")
-  MESSAGE(STATUS "  NetCDF_INC:       " "${NetCDF_INCLUDE_DIR}")
-  MESSAGE(STATUS "  NetCDF_LIBS:     " "${NetCDF_LIBRARIES}")
+  MESSAGE(STATUS "  netCDF:           " "${NETCDF_FOUND}")
+  MESSAGE(STATUS "  netCDF_INC:       " "${NETCDF_INCLUDE_DIR}")
+  MESSAGE(STATUS "  netCDF_LIBS:      " "${NETCDF_LIBRARIES}")
 ELSE()
-  MESSAGE(STATUS "Library not found: >NetCDF_FOUND< ")
-  MESSAGE(WARNING "  \n Missing: >NetCDF_INCLUDE_DIR< , >NetCDF_LIBRARY<,  >NetCDFF_LIBRARY< \n some functionalities will be disabled")
-ENDIF(NetCDF_FOUND)
+  MESSAGE(STATUS "Library not found: netCDF ")
+  MESSAGE(WARNING "  \n Missing: <NETCDF_INCLUDE_DIR> , <NETCDF_LIBRARY>,  <NETCDFF_LIBRARY> \n some functionalities will be disabled")
+ENDIF()
 
 # ---------------------- #
 # -- HDF5 libraries -- #
@@ -60,39 +60,39 @@ ENDIF()
 SET(CMAKE_Fortran_MODULE_DIRECTORY
   ${PROJECT_BINARY_DIR}/fmodules CACHE PATH "Directory for Fortran modules")
 
-SET(ElmerIce_SRC ElmerIceUtils.F90 AIFlowSolve_nlD2.F90 AIFlowSolve_nlS2.F90 
-    CaffeSolver.F90 ComputeDevStress.F90 ComputeEigenValues.F90 
-    ComputeNormal.F90 ComputeStrainRate.F90 DeformationalHeat.F90 
-    EPLSolver.F90 FabricSolve.F90 Flowdepth.F90 
-    ForceToStress.F90 GetHydrostaticLoads.F90  GolfLaw.F90 
-    GroundedSolver.F90  IntegratedVelocity.F90  IDSSolver.F90 
-    PorousSolve.F90  pointwise.F90  SIASolver.F90  SSASolver.F90  
-    ThicknessSolver.F90  TemperateIce.F90  ExportVertically.F90  
-    AdjointSolver.F90  DJDBeta_Adjoint.F90  DJDmu_Adjoint.F90  
-    CostSolver_Adjoint.F90  DJDBeta_Robin.F90  DJDmu_Robin.F90  
-    CostSolver_Robin.F90  m1qn3.F  Grid2DInterpolator.F90  
+SET(ElmerIce_SRC ElmerIceUtils.F90 AIFlowSolve_nlD2.F90 AIFlowSolve_nlS2.F90
+    CaffeSolver.F90 ComputeDevStress.F90 ComputeEigenValues.F90
+    ComputeNormal.F90 ComputeStrainRate.F90 DeformationalHeat.F90
+    EPLSolver.F90 FabricSolve.F90 Flowdepth.F90
+    ForceToStress.F90 GetHydrostaticLoads.F90  GolfLaw.F90
+    GroundedSolver.F90  IntegratedVelocity.F90  IDSSolver.F90
+    PorousSolve.F90  pointwise.F90  SIASolver.F90  SSASolver.F90
+    ThicknessSolver.F90  TemperateIce.F90  ExportVertically.F90
+    AdjointSolver.F90  DJDBeta_Adjoint.F90  DJDmu_Adjoint.F90
+    CostSolver_Adjoint.F90  DJDBeta_Robin.F90  DJDmu_Robin.F90
+    CostSolver_Robin.F90  m1qn3.F  Grid2DInterpolator.F90
     Optimize_m1qn3Parallel.F90  OutputStrainHeating.F90  UpdateExport.F90
     IntegrateVertically.F90 EnthalpySolver.F90 SubShelfMelt.F90
     ./Adjoint/Adjoint_LinearSolver.F90 ./Adjoint/Adjoint_CostDiscSolver.F90
-    ./Adjoint/Adjoint_CostContSolver.F90 ./Adjoint/Adjoint_CostRegSolver.F90 
+    ./Adjoint/Adjoint_CostContSolver.F90 ./Adjoint/Adjoint_CostRegSolver.F90
     ./Adjoint/Adjoint_GradientValidation.F90
-    ./AdjointStokes/AdjointStokes_GradientMu.F90 
-    ./AdjointStokes/AdjointStokes_GradientBetaSolver.F90  
+    ./AdjointStokes/AdjointStokes_GradientMu.F90
+    ./AdjointStokes/AdjointStokes_GradientBetaSolver.F90
     ./AdjointSSA/AdjointSSA_AdjointSolver.F90 ./AdjointSSA/AdjointSSA_CostDiscSolver.F90
     ./AdjointSSA/AdjointSSA_CostRegSolver.F90 ./AdjointSSA/AdjointSSA_SSASolver.F90
     ./AdjointSSA/AdjointSSA_CostContSolver.F90 ./AdjointSSA/AdjointSSA_CostFluxDivSolver.F90
     ./AdjointSSA/AdjointSSA_CostTaubSolver.F90
-    ./AdjointSSA/AdjointSSA_GradientSolver.F90 
+    ./AdjointSSA/AdjointSSA_GradientSolver.F90
     ./AdjointThickness/AdjointThickness_GradientSolver.F90   ./AdjointThickness/AdjointThickness_ThicknessSolver.F90
     ./Permafrost/PermafrostMaterials.F90 ./Permafrost/Permafrost_Utils.F90 ./Permafrost/Permafrost_HTEQ.F90
     ./Permafrost/Permafrost_Darcy.F90 ./Permafrost/Permafrost_solute.F90 ./Permafrost/Permafrost_solid.F90
     SurfaceBoundaryEnthalpy.F90
-    Calving.F90 FrontDisplacement.F90 
-    TwoMeshes.F90 ProjectCalving.F90 ComputeCalvingNormal.F90 
+    Calving.F90 FrontDisplacement.F90
+    TwoMeshes.F90 ProjectCalving.F90 ComputeCalvingNormal.F90
     CalvingGeometry.F90 Calving3D.F90 Calving3D_lset.F90 CalvingGlacierAdvance3D.F90 CalvingRemesh.F90
-    CalvingFrontAdvance3D.F90 Emergence.F90 SSAmask.F90 
+    CalvingFrontAdvance3D.F90 Emergence.F90 SSAmask.F90
     GlaDSCoupledSolver.F90 GlaDSchannelSolver.F90 Flotation.F90
-    BasalMelt3D.F90 CalvingHydroInterp.F90 HydroRestart.F90 
+    BasalMelt3D.F90 CalvingHydroInterp.F90 HydroRestart.F90
     GMValid.F90 Scalar_OUTPUT_Glacier.F90 IcyMaskSolver.F90
     Weertman2Coulomb.F90)
 
@@ -123,7 +123,7 @@ ENDIF()
 ADD_LIBRARY(ElmerIceSolvers SHARED ${ElmerIce_SRC})
 
 # Library object
-SET_TARGET_PROPERTIES(ElmerIceSolvers PROPERTIES PREFIX "")  
+SET_TARGET_PROPERTIES(ElmerIceSolvers PROPERTIES PREFIX "")
 SET_TARGET_PROPERTIES(ElmerIceSolvers PROPERTIES
     LINKER_LANGUAGE Fortran
     LIBRARY_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/fem/src/modules
@@ -137,7 +137,7 @@ ENDIF()
 TARGET_LINK_LIBRARIES(ElmerIceSolvers Elmer::MPI_Fortran elmersolver ElmerIceUtils)
 
 IF(HAVE_NETCDF)
-TARGET_LINK_LIBRARIES(ElmerIceSolvers ${NETCDF_LIBRARIES})
+  TARGET_LINK_LIBRARIES(ElmerIceSolvers ${NETCDF_LIBRARIES})
 ENDIF()
 IF(HAVE_HDF5)
 TARGET_LINK_LIBRARIES(ElmerIceSolvers ${PHDF5_LIBRARIES})
@@ -164,6 +164,6 @@ INSTALL(TARGETS ElmerIceSolvers LIBRARY DESTINATION "share/elmersolver/lib"
 IF(HAVE_NETCDF)
   ADD_SUBDIRECTORY(GridDataReader)
 ENDIF()
-IF(WITH_ScatteredDataInterpolator)	
+IF(WITH_ScatteredDataInterpolator)
   ADD_SUBDIRECTORY(ScatteredDataInterpolator)
 ENDIF()

elmerice/Tests/FrictionHeatMasked/DummySolver.f90

@@ -31,6 +31,7 @@ RECURSIVE SUBROUTINE DummySolver( Model,Solver,Timestep,TransientSimulation )
      VarPerm => Var % Perm
      VarDOFs =  Var % DOFs
      VarValues => Var % Values
+     VarValues = 0.0_dp
   ELSE
      CALL FATAL('DummySolver','No Variable associated')
   END IF
@@ -41,13 +42,13 @@ RECURSIVE SUBROUTINE DummySolver( Model,Solver,Timestep,TransientSimulation )
  !       VarValues(VarDOFs*(VarPerm(i) - 1)+j) = k
  !    END DO
  ! END DO
- ! VarValues = 0.0_dp
   DO t=1, Solver % Mesh % NumberOfBoundaryElements
       ! get element information
       Element => GetBoundaryElement(t)
       IF ( .NOT.ActiveBoundaryElement() ) CYCLE
       BC => GetBC()
       n = GetElementNOFNodes()
+      IF(ANY(VarPerm(Element % NodeIndexes) == 0)) CYCLE
       VarValues(VarPerm(Element % NodeIndexes)) = ListGetReal(BC,TRIM(Solver % Variable % Name),n,Element % NodeIndexes,GotIt)
    END DO
 

elmerice/Tests/FrictionHeatMasked/test.sif

@@ -180,6 +180,10 @@ Solver 6
   Procedure = File "DummySolver" "DummySolver"
   Variable =  -dofs 1 "Friction Load"
   Exported Variable 1 = -dofs 1 "Friction Load Mask"
+
+! This seems to result to NaN for some reason on some platforms.
+! So supress it for now.
+  Exec Solver = never
 End
 
 !creates output