fix Main_DOODZ

tduretz · Feb 19, 2024 · 3649b7a · 3649b7a
1 parent 868ee8b
commit 3649b7a
Show file tree

Hide file tree

Showing 2 changed files with 8 additions and 6 deletions.
diff --git a/MDLIB/Main_DOODZ.c b/MDLIB/Main_DOODZ.c
@@ -111,8 +111,7 @@ void RunMDOODZ(char *inputFileName, MdoodzSetup *setup) {
     markers particles = PartAlloc(inputFile.particles, &input.model );
 
     // Allocate marker chain
-    markers       topo_chain, topo_chain_ini;
-    topo_chain.Nb_part = 0;
+    markers      topo_chain, topo_chain_ini;
     surface      topo, topo_ini;
     if (input.model.free_surface == 1 ) AllocateMarkerChain( &topo,     &topo_chain, input.model );
     if (input.model.free_surface == 1 ) AllocateMarkerChain( &topo_ini, &topo_chain_ini, input.model );
@@ -137,7 +136,8 @@ void RunMDOODZ(char *inputFileName, MdoodzSetup *setup) {
         }
 
         // Initial grid tags
-        SetBCs(*setup->SetBCs, &input, &mesh, &topo_chain);
+        SetBCs(*setup->SetBCs, &input, &mesh, &topo);
+
         if (input.model.free_surface == 1 ) {
 
             // Define the horizontal position of the surface marker chain
@@ -199,7 +199,7 @@ void RunMDOODZ(char *inputFileName, MdoodzSetup *setup) {
           printf("Running with crazy conductivity for the asthenosphere!!\n");
         }
         // Initial solution fields
-        SetBCs(*setup->SetBCs, &input, &mesh, &topo_chain);
+        SetBCs(*setup->SetBCs, &input, &mesh, &topo);
 
         if (input.model.mechanical == 1) {
             InitialiseSolutionFields( &mesh, &input.model );
@@ -225,7 +225,7 @@ void RunMDOODZ(char *inputFileName, MdoodzSetup *setup) {
         P2Mastah( &input.model, particles, input.materials.Cv,    &mesh, mesh.Cv, mesh.BCp.type,  0, 0, interp, cent, 1);
         P2Mastah( &input.model, particles, input.materials.Qr,    &mesh, mesh.Qr, mesh.BCp.type,  0, 0, interp, cent, 1);
 
-        SetBCs(*setup->SetBCs, &input, &mesh, &topo_chain);
+        SetBCs(*setup->SetBCs, &input, &mesh, &topo);
         if (input.model.initial_cooling == 1 ) ThermalSteps( &mesh, input.model,  mesh.rhs_t, &particles, input.model.cooling_duration, input.scaling );
         if (input.model.therm_perturb == 1 ) SetThermalPert( &mesh, input.model, input.scaling );
         Interp_Grid2P_centroids2( particles, particles.T,    &mesh, mesh.T, mesh.xvz_coord,  mesh.zvx_coord,  mesh.Nx-1, mesh.Nz-1, mesh.BCt.type, &input.model );
@@ -593,7 +593,7 @@ void RunMDOODZ(char *inputFileName, MdoodzSetup *setup) {
                 printf("Running with normal conductivity for the asthenosphere...\n");
             }
             // Allocate and initialise solution and RHS vectors
-            SetBCs(*setup->SetBCs, &input, &mesh, &topo_chain);
+            SetBCs(*setup->SetBCs, &input, &mesh, &topo);
 
             // Reset fields and BC values if needed
             //        if ( input.model.pure_shear_ALE == 1 ) InitialiseSolutionFields( &mesh, &input.model );

diff --git a/SETS/RiverTomCrouteinf.c b/SETS/RiverTomCrouteinf.c
@@ -131,6 +131,8 @@ SetBC SetBCVz(MdoodzInput *instance, POSITION position, Coordinates coordinates)
   const double Vz_corr   = surf_Winc*surf_Vinc / Lx;
   const double VzS = (-0.5 * V_tot * (instance->topo_height->west - instance->model.zmin) / Lx) -(0.5 * V_tot * (instance->topo_height->east - instance->model.zmin) / Lx) + Vz_corr;
 
+// printf("topo_height->west = %2.2e --- topo_height->east = %2.2e ", instance->topo_height->west, instance->topo_height->east);
+
   // Set boundary nodes types and values
   if (position == W || position == SW || position == NW ) {
     bc.type  = constant_shear_stress;