From 18f373a193af3089c490c7a74edba21293b8969f Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Mon, 28 May 2018 16:45:55 +0200 Subject: [PATCH 1/9] Removed three unused global functions --- GEO/geo_sfc.cpp | 242 ------------------------------------------------ GEO/geo_sfc.h | 4 +- 2 files changed, 1 insertion(+), 245 deletions(-) diff --git a/GEO/geo_sfc.cpp b/GEO/geo_sfc.cpp index 055509e35..a5cb9dc3f 100644 --- a/GEO/geo_sfc.cpp +++ b/GEO/geo_sfc.cpp @@ -892,66 +892,7 @@ void Surface::CreateTIN(void) if (TIN->Triangles.size() == 0) delete TIN; } -/************************************************************************** - GeoLib-Method: - Task: - Programing: - 03/2004 OK Implementation - last modification:08/2005 CC -**************************************************************************/ -void GEOWriteSurfaceTINs(string file_path) -{ - vector::iterator ps = surface_vector.begin(); - Surface* m_surface = NULL; - while (ps != surface_vector.end()) - { - m_surface = *ps; - if (!m_surface->TIN) - { - ++ps; - continue; - } - else - { - m_surface->WriteTIN(file_path); - ++ps; - } - } -} -/************************************************************************** - GeoLib-Method: - Task: - Programing: - 03/2004 OK Implementation - last modification: -**************************************************************************/ -void GEOWriteSurfaceTINsTecplot(string file_path) -{ - vector::iterator ps = surface_vector.begin(); - Surface* m_surface = NULL; - while (ps != surface_vector.end()) - { - m_surface = *ps; - if (!m_surface->TIN) - { - ++ps; - continue; - } - else - { - m_surface->WriteTINTecplot(file_path); - ++ps; - } - } -} - -/* - void Surface::operator = (const Surface m_surface) - { - name = m_surface.name; - } - */ /************************************************************************** GEOLib-Method: @@ -1256,71 +1197,6 @@ void Surface::WriteTIN(const std::string& file_path) // tin_file.close(); } -/************************************************************************** - GeoLib-Method: - Task: - Programing: - 03/2004 OK Implementation - last modification:08/2005 CC -**************************************************************************/ -void GEOCreateLayerSurfaceTINs(int nb_prism_layers, double thickness_prism_layers) -{ - long i, j; - long no_TIN_Triangles; - vector::iterator ps = surface_vector.begin(); - Surface* m_surface = NULL; - Surface* m_surface_new = NULL; - char layer_number[10]; - CTriangle* m_triangle = NULL; - vector new_surfaces; - - // Go through the existing surfaces (GEO_type) - while (ps != surface_vector.end()) - { - m_surface = *ps; - if (!m_surface->TIN) - { - ++ps; - continue; - } - no_TIN_Triangles = (long)m_surface->TIN->Triangles.size(); - if (m_surface->type == 0) - // Go through the layers - for (i = 0; i < nb_prism_layers; i++) - { - m_surface_new = new Surface; - sprintf(layer_number, "%ld", i + 1); - m_surface_new->name = m_surface->name + "_L" + layer_number; - m_surface_new->type = m_surface->type; - m_surface_new->TIN = new CTIN; - m_surface_new->TIN->name = m_surface_new->name; - for (j = 0; j < no_TIN_Triangles; j++) - { - m_triangle = new CTriangle; - m_triangle->number = m_surface->TIN->Triangles[j]->number; - m_triangle->x[0] = m_surface->TIN->Triangles[j]->x[0]; - m_triangle->x[1] = m_surface->TIN->Triangles[j]->x[1]; - m_triangle->x[2] = m_surface->TIN->Triangles[j]->x[2]; - m_triangle->y[0] = m_surface->TIN->Triangles[j]->y[0]; - m_triangle->y[1] = m_surface->TIN->Triangles[j]->y[1]; - m_triangle->y[2] = m_surface->TIN->Triangles[j]->y[2]; - m_triangle->z[0] = m_surface->TIN->Triangles[j]->z[0] - (i + 1) * thickness_prism_layers; - m_triangle->z[1] = m_surface->TIN->Triangles[j]->z[1] - (i + 1) * thickness_prism_layers; - m_triangle->z[2] = m_surface->TIN->Triangles[j]->z[2] - (i + 1) * thickness_prism_layers; - m_surface_new->TIN->Triangles.push_back(m_triangle); - } - new_surfaces.push_back(m_surface_new); - } - ++ps; - } - // - int new_surfaces_vector_length = (int)new_surfaces.size(); - for (i = 0; i < new_surfaces_vector_length; i++) - { - m_surface_new = new_surfaces[i]; - surface_vector.push_back(m_surface_new); - } -} /************************************************************************** GeoLib-Method: Task: @@ -1372,124 +1248,6 @@ void Surface::WriteTINTecplot(const std::string& file_path) } } -/************************************************************************** - GeoLib-Method: - Task: - Programing: - 04/2004 OK Implementation - 01/2005 OK TIN case - last modification: -**************************************************************************/ -/*vector Surface::GetMSHNodesClose() - { - long i,j; - long no_nodes = 0; - long *nodes_array = NULL; - long no_points; - CGLPoint m_node; - CGLPolyline* m_polyline = NULL; - CGLPolyline* m_polyline1 = NULL; - CGLPolyline* m_polyline2 = NULL; - vectornodes_vector; - double xp[3],yp[3],zp[3]; - CTriangle *m_triangle = NULL; - long no_triangles; - - list::const_iterator p = polyline_of_surface_list.begin(); - nodes_vector.clear(); - - switch(data_type) { - //-------------------------------------------------------------------- - case 0: // surface polygon - nodes_array = GetPointsIn(this,&no_nodes);//change this into m_sfc todo - for(i=0;iTriangles.size(); - for(i=0;iTriangles[i]; - xp[0] = m_triangle->x[0]; - yp[0] = m_triangle->y[0]; - zp[0] = m_triangle->z[0]; - xp[1] = m_triangle->x[1]; - yp[1] = m_triangle->y[1]; - zp[1] = m_triangle->z[1]; - xp[2] = m_triangle->x[2]; - yp[2] = m_triangle->y[2]; - zp[2] = m_triangle->z[2]; - for(j=0;jpoint_vector.size(); - - for(j=0;jpoint_vector[i]->x; - yp[0] = m_polyline1->point_vector[i]->y; - zp[0] = m_polyline1->point_vector[i]->z; - xp[1] = m_polyline1->point_vector[i+1]->x; - yp[1] = m_polyline1->point_vector[i+1]->y; - zp[1] = m_polyline1->point_vector[i+1]->z; - xp[2] = m_polyline2->point_vector[i]->x; - yp[2] = m_polyline2->point_vector[i]->y; - zp[2] = m_polyline2->point_vector[i]->z; - if(m_node.IsInsideTriangle(xp,yp,zp)) { - nodes_vector.push_back(nodes_array[j]); - } - // second triangle of quad - xp[0] = m_polyline2->point_vector[i]->x; - yp[0] = m_polyline2->point_vector[i]->y; - zp[0] = m_polyline2->point_vector[i]->z; - xp[1] = m_polyline2->point_vector[i+1]->x; - yp[1] = m_polyline2->point_vector[i+1]->y; - zp[1] = m_polyline2->point_vector[i+1]->z; - xp[2] = m_polyline1->point_vector[i+1]->x; - yp[2] = m_polyline1->point_vector[i+1]->y; - zp[2] = m_polyline1->point_vector[i+1]->z; - if(m_node.IsInsideTriangle(xp,yp,zp)) { - nodes_vector.push_back(nodes_array[j]); - } - } // no_points - } // no_nodes - - break; - } - return nodes_vector; - }*/ - /************************************************************************** GeoLib-Method: Task: diff --git a/GEO/geo_sfc.h b/GEO/geo_sfc.h index 679426ee3..6c3ac3b14 100644 --- a/GEO/geo_sfc.h +++ b/GEO/geo_sfc.h @@ -120,9 +120,7 @@ extern void GEOSurfaceTopology(void); extern void GEOUnselectSFC(); // OK // TIN #define TIN_FILE_EXTENSION ".tin" -extern void GEOWriteSurfaceTINs(const std::string&); // TIN -extern void GEOCreateLayerSurfaceTINs(int, double); // TIN -extern void GEOWriteSurfaceTINsTecplot(const std::string&); + extern int sfc_ID_max; // MSH void MSHUnselectSFC(); // OK From 4909df239b5efb048a08317400063402d22a8782 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Mon, 28 May 2018 16:47:12 +0200 Subject: [PATCH 2/9] Added the forgot "break" in some switch --- FEM/rf_mfp_new.cpp | 1 + FEM/rf_mmp_new.cpp | 1 + FEM/rf_msp_new.cpp | 1 + FEM/rf_pcs.cpp | 5 +++-- FEM/rfmat_cp.cpp | 2 ++ 5 files changed, 8 insertions(+), 2 deletions(-) diff --git a/FEM/rf_mfp_new.cpp b/FEM/rf_mfp_new.cpp index f43ebb62a..5dabf02aa 100644 --- a/FEM/rf_mfp_new.cpp +++ b/FEM/rf_mfp_new.cpp @@ -1211,6 +1211,7 @@ double CFluidProperties::Density(double* variables) std::cout << " Error - Density Model 18 not implemented, usind dummy density of 1000." << "\n"; density = 1000; // Achtung - dummy + break; case 19: // KG44 get the density from GEMS calculations // seems complicated, as we probably have to call GEMS.....or take values from last GEMS calculation diff --git a/FEM/rf_mmp_new.cpp b/FEM/rf_mmp_new.cpp index 0e3774cdc..9b5928bc6 100644 --- a/FEM/rf_mmp_new.cpp +++ b/FEM/rf_mmp_new.cpp @@ -978,6 +978,7 @@ std::ios::pos_type CMediumProperties::Read(std::ifstream* mmp_file) // if -1, threshold is constant in >> permeability_strain_model_value[4]; // lower limit in >> permeability_strain_model_value[5]; // uper limit + break; default: cout << "Error in MMPRead: no valid permeability strain model" << "\n"; diff --git a/FEM/rf_msp_new.cpp b/FEM/rf_msp_new.cpp index e7e000dc3..1b78c4430 100644 --- a/FEM/rf_msp_new.cpp +++ b/FEM/rf_msp_new.cpp @@ -8168,6 +8168,7 @@ double CSolidProperties::E_Function(int dim, const ElementValue_DM* ele_val, int } CalPrinStrDir(stress, prin_str, prin_dir, size); return_value = GetCurveValue((int)E_Function_Model_Value[0], 0, prin_str[0], &valid); + break; } default: return_value = 1.; diff --git a/FEM/rf_pcs.cpp b/FEM/rf_pcs.cpp index f602fc1fc..6c3920a7a 100644 --- a/FEM/rf_pcs.cpp +++ b/FEM/rf_pcs.cpp @@ -12615,10 +12615,11 @@ void CRFProcess::AssembleParabolicEquationRHSVector(CNode* m_nod) */ } } - if (m_edg->GetMark()) - break; //---------------------------------------------------------------- // ToDo + //if (m_edg->GetMark()) + // break; + break; default: cout << "Warning in CRFProcess::AssembleParabolicEquationRHSVector - not implemented for this element " "type" diff --git a/FEM/rfmat_cp.cpp b/FEM/rfmat_cp.cpp index 9ad983aeb..5fb41517c 100644 --- a/FEM/rfmat_cp.cpp +++ b/FEM/rfmat_cp.cpp @@ -1261,6 +1261,7 @@ int CompProperties::GetNumberDiffusionValuesCompProperties(int diffusion_model) break; /* Keine Diffusion */ case 0: /* curve */ n = 1; + break; case 1: n = 1; break; /* Konstanter Diffusionswert */ @@ -1379,6 +1380,7 @@ int CompProperties::GetNumberIsothermValuesCompProperties(int isotherm) { case -1: n = 0; /* no isotherm */ + break; case 0: n = 1; break; /* get KD from curve derivative */ From 132d7c68b177ce0dbb18e1eeff24e3a6c6563ea8 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Tue, 29 May 2018 10:27:25 +0200 Subject: [PATCH 3/9] Fixed some compilation warning --- FEM/Output.cpp | 9 --------- FEM/fem_ele_vec.cpp | 28 +++++++--------------------- FEM/rf_pcs.cpp | 15 ++++----------- FEM/rf_random_walk.cpp | 1 - FEM/rf_react.cpp | 4 ++-- 5 files changed, 13 insertions(+), 44 deletions(-) diff --git a/FEM/Output.cpp b/FEM/Output.cpp index 98897b392..397be5e66 100644 --- a/FEM/Output.cpp +++ b/FEM/Output.cpp @@ -3477,7 +3477,6 @@ void COutput::PCONWriteDOMDataTEC() void COutput::WriteTECNodePCONData(fstream& tec_file) { const size_t nName(_pcon_value_vector.size()); - int nidx_dm[3]; std::vector PconIndex(nName); // m_msh = GetMSH(); @@ -3509,14 +3508,6 @@ void COutput::WriteTECNodePCONData(fstream& tec_file) // XYZ double x[3] = {m_msh->nod_vector[j]->getData()[0], m_msh->nod_vector[j]->getData()[1], m_msh->nod_vector[j]->getData()[2]}; - // x[0] = m_msh->nod_vector[j]->X(); - // x[1] = m_msh->nod_vector[j]->Y(); - // x[2] = m_msh->nod_vector[j]->Z(); - // Amplifying DISPLACEMENTs - if (M_Process || MH_Process) // WW - - for (size_t k = 0; k < max_dim + 1; k++) - x[k] += out_amplifier * m_pcs->GetNodeValue(m_msh->nod_vector[j]->GetIndex(), nidx_dm[k]); for (size_t i = 0; i < 3; i++) tec_file << x[i] << " "; // NOD values diff --git a/FEM/fem_ele_vec.cpp b/FEM/fem_ele_vec.cpp index b2116393a..a66ec10b8 100644 --- a/FEM/fem_ele_vec.cpp +++ b/FEM/fem_ele_vec.cpp @@ -3790,7 +3790,13 @@ void CFiniteElementVec::GlobalAssembly_RHS() ---------------------- -----------------------------------------------------------------*/ ElementValue_DM::ElementValue_DM(CElem * ele, const int NGP, bool HM_Staggered) - : NodesOnPath(NULL), orientation(NULL), Strain(NULL), Strain_t_ip(NULL) + : Stress(NULL), Stress_i(NULL), Stress_j(NULL), pStrain(NULL), + y_surface(NULL), prep0(NULL), e_i(NULL), xi(NULL), MatP(NULL), + Strain_Kel(NULL), Strain_Max(NULL), Strain_pl(NULL), + Strain_t_ip(NULL), e_pl(NULL), ev_loc_nr_res(NULL), + lambda_pl(NULL), Strain(NULL), NodesOnPath(NULL), + orientation(NULL), scalar_aniso_comp(NULL), + scalar_aniso_tens(NULL) { int Plastic = 1; const int LengthMat = 7; // Number of material parameter of SYS model. @@ -3799,24 +3805,10 @@ void CFiniteElementVec::GlobalAssembly_RHS() CSolidProperties* sdp = NULL; int ele_dim; // - Stress = NULL; - pStrain = NULL; - prep0 = NULL; - e_i = NULL; - xi = NULL; - MatP = NULL; - NodesOnPath = NULL; - orientation = NULL; MshElemType::type ele_type = ele->GetElementType(); ele_dim = ele->GetDimension(); sdp = msp_vector[ele->GetPatchIndex()]; Plastic = sdp->Plastictity(); - Strain_Kel = NULL; - Strain_Max = NULL; - Strain_pl = NULL; - e_pl = NULL; - lambda_pl = NULL; - ev_loc_nr_res = NULL; if (ele_dim == 2) LengthBS = 4; @@ -3837,8 +3829,6 @@ void CFiniteElementVec::GlobalAssembly_RHS() Stress = Stress_i; if (HM_Staggered) Stress_j = new Matrix(LengthBS, NGPoints); - else - Stress_j = NULL; // for HM coupling iteration // if (Plastic > 0) { @@ -3847,8 +3837,6 @@ void CFiniteElementVec::GlobalAssembly_RHS() *y_surface = 0.0; *pStrain = 0.0; } - else - y_surface = NULL; *Stress = 0.0; if (Plastic == 2) // Rotational hardening model @@ -3913,8 +3901,6 @@ void CFiniteElementVec::GlobalAssembly_RHS() disp_j = 0.0; tract_j = 0.0; Localized = false; - scalar_aniso_comp = NULL; // WX: 11.2011 plasticity bedding - scalar_aniso_tens = NULL; if (sdp->Plasticity_Bedding) { scalar_aniso_comp = new Matrix(NGPoints); diff --git a/FEM/rf_pcs.cpp b/FEM/rf_pcs.cpp index 6c3920a7a..3db1fd9db 100644 --- a/FEM/rf_pcs.cpp +++ b/FEM/rf_pcs.cpp @@ -8647,21 +8647,14 @@ void CRFProcess::CalcSecondaryVariables(bool initial) case FiniteElement::TWO_PHASE_FLOW: break; case FiniteElement::RICHARDS_FLOW: // Richards flow - // WW + case FiniteElement::MULTI_PHASE_FLOW: + case FiniteElement::DEFORMATION_H2: // H2M CalcSecondaryVariablesUnsaturatedFlow(initial); - break; - case FiniteElement::DEFORMATION || FiniteElement::DEFORMATION_FLOW || FiniteElement::DEFORMATION_DYNAMIC: - if (type == 42) // H2M //WW - CalcSecondaryVariablesUnsaturatedFlow(initial); - break; case FiniteElement::PS_GLOBAL: CalcSecondaryVariablesPSGLOBAL(); // WW break; default: - if (type == 1212) - // WW - CalcSecondaryVariablesUnsaturatedFlow(initial); break; } } @@ -11834,7 +11827,7 @@ void CRFProcess::CalcELEFluxes(const GEOLIB::Polyline* const ply, double* result for (int k = 0; k < int(vec_nodes_edge.size()); k++) { Point_on_Geo = false; - for (int l = 0; l < int(nod_vector_at_geo.size()); l++) + for (std::size_t l = 0; l < nod_vector_at_geo.size(); l++) { if (vec_nodes_edge[k]->GetIndex() == nod_vector_at_geo[l]) { @@ -12123,7 +12116,7 @@ void CRFProcess::CalcELEMassFluxes(const GEOLIB::Polyline* const ply, std::strin for (int k = 0; k < int(vec_nodes_edge.size()); k++) { Point_on_Geo = false; - for (int l = 0; l < int(nod_vector_at_geo.size()); l++) + for (std::size_t l = 0; l < nod_vector_at_geo.size(); l++) { if (vec_nodes_edge[k]->GetIndex() == nod_vector_at_geo[l]) { diff --git a/FEM/rf_random_walk.cpp b/FEM/rf_random_walk.cpp index 4e7f03746..dcfb289a8 100644 --- a/FEM/rf_random_walk.cpp +++ b/FEM/rf_random_walk.cpp @@ -6133,7 +6133,6 @@ void DATWriteParticleControlPlaneFile(int current_time_step, string control_plan if (pos==0) vtk_file << "id, starting_time, arrival_time_step_time, arrival-time-cp , start_coor_x, start_coor_y, start_coor_z, control_plane_x, control_plane_y, control_plane_z, time_correction_factor, x_back_shift " << endl; - int id = 0; double starting_time = 0.0; double time_correction = 0.0; double arrival_time_step_time = 0.0; diff --git a/FEM/rf_react.cpp b/FEM/rf_react.cpp index fea24b915..776b70112 100644 --- a/FEM/rf_react.cpp +++ b/FEM/rf_react.cpp @@ -427,12 +427,12 @@ void REACT::ExecutePQCString(void) ii = 0; int idx, idy; bool firstinput; - char string[4]; //, string1[4], string2[4]; +// char string[4]; //, string1[4], string2[4]; for (std::size_t j = 0; j < rankrankliststore[myrank].size(); j++) { firstinput = true; idy = rankrankliststore[myrank][j]; - sprintf(string, "%li", static_cast(idy)); + //sprintf(string, "%li", static_cast(idy)); for (std::size_t i = 0; i < ranknodeliststore[idy].size(); i++) { idx = ranknodeliststore[idy][i]; From deeb218e567ea3d32f4f9be1b90b4eef6b4c80ac Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Tue, 29 May 2018 11:06:42 +0200 Subject: [PATCH 4/9] Fixed some compilation warnings of "ignoring return value" --- Base/FileTools.cpp | 3 ++- FEM/fem_ele_std.cpp | 2 +- FEM/files0.cpp | 3 ++- FEM/problem.cpp | 6 +++--- FEM/rf_mfp_new.cpp | 2 +- FEM/rf_pcs.cpp | 6 +++--- 6 files changed, 12 insertions(+), 10 deletions(-) diff --git a/Base/FileTools.cpp b/Base/FileTools.cpp index 7fe74b9b5..9add2f9b0 100644 --- a/Base/FileTools.cpp +++ b/Base/FileTools.cpp @@ -155,7 +155,8 @@ std::string getCwd() #ifdef WIN32 _getcwd(cwd, FILENAME_MAX); #else - getcwd(cwd, FILENAME_MAX); + char* unused __attribute__((unused)); + unused = getcwd(cwd, FILENAME_MAX); #endif return cwd; diff --git a/FEM/fem_ele_std.cpp b/FEM/fem_ele_std.cpp index c4b840558..706a46299 100644 --- a/FEM/fem_ele_std.cpp +++ b/FEM/fem_ele_std.cpp @@ -6778,7 +6778,7 @@ string CFiniteElementStd::Cal_GP_Velocity_DuMux(int* i_ind, CRFProcess* m_pcs, i { cout << "The program is canceled because there is a phase used which is not considered yet!" << "\n"; - system("Pause"); + //system("Pause"); exit(0); } } diff --git a/FEM/files0.cpp b/FEM/files0.cpp index 9cad05fa3..c9947cee0 100644 --- a/FEM/files0.cpp +++ b/FEM/files0.cpp @@ -718,7 +718,8 @@ char* ReadString(void) { char* s = (char*)malloc(256); // char *s = new char[256];//CC - scanf(" %s%*[^\n]%*c", s); + int unused __attribute__((unused)); + unused = scanf(" %s%*[^\n]%*c", s); // int a = (int)strlen(s); // delete[] s; // s = new char[a+1];//CC diff --git a/FEM/problem.cpp b/FEM/problem.cpp index 5c7f88d54..e7bcada51 100644 --- a/FEM/problem.cpp +++ b/FEM/problem.cpp @@ -1982,7 +1982,7 @@ inline double Problem::MultiPhaseFlow() { std::cout << "Error running Eclipse!" << "\n"; - system("Pause"); + //system("Pause"); exit(0); } } @@ -2016,7 +2016,7 @@ inline double Problem::MultiPhaseFlow() << "\n"; std::cout << "The run is terminated now ..." << "\n"; - system("Pause"); + //system("Pause"); exit(0); } FluidProp = MFPGet("GAS"); @@ -2027,7 +2027,7 @@ inline double Problem::MultiPhaseFlow() << "\n"; std::cout << "The run is terminated now ..." << "\n"; - system("Pause"); + //system("Pause"); exit(0); } } diff --git a/FEM/rf_mfp_new.cpp b/FEM/rf_mfp_new.cpp index 5dabf02aa..a3a71fbdd 100644 --- a/FEM/rf_mfp_new.cpp +++ b/FEM/rf_mfp_new.cpp @@ -1318,7 +1318,7 @@ double CFluidProperties::GetElementValueFromNodes(long ElementIndex, int GPIndex << "\n"; cout << "The run is terminated now ..." << "\n"; - system("Pause"); + //system("Pause"); exit(0); } diff --git a/FEM/rf_pcs.cpp b/FEM/rf_pcs.cpp index 3db1fd9db..497fe5866 100644 --- a/FEM/rf_pcs.cpp +++ b/FEM/rf_pcs.cpp @@ -15270,7 +15270,7 @@ void CRFProcess::CalculateFluidDensitiesAndViscositiesAtNodes(CRFProcess* m_pcs) { cout << "Density calculation of water was not possible" << "\n"; - system("Pause"); + //system("Pause"); exit(0); } // calculate new moles of H2O @@ -15353,7 +15353,7 @@ void CRFProcess::CalculateFluidDensitiesAndViscositiesAtNodes(CRFProcess* m_pcs) { cout << "Density calculation of gas was not possible" << "\n"; - system("Pause"); + //system("Pause"); exit(0); } // set new density to nodes @@ -15695,7 +15695,7 @@ void CRFProcess::Phase_Transition_CO2(CRFProcess* m_pcs, int Step) cout << "The volume is not equal before and after the calculation of CO2 phase transition! " << "\n"; cout << "Before: " << Volume_eff << " After: " << gas.volume + liquid.volume << "\n"; - system("Pause"); + //system("Pause"); exit(0); } From ee231f6ea87d8444bca20d165f770c59b093d240 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Tue, 29 May 2018 13:12:05 +0200 Subject: [PATCH 5/9] Removed unused source code in rf.cpp --- Base/memory.cpp | 182 ------------------------------------------------ Base/memory.h | 13 ---- OGS/rf.cpp | 28 -------- 3 files changed, 223 deletions(-) diff --git a/Base/memory.cpp b/Base/memory.cpp index b6c4b0f3b..160c077a9 100644 --- a/Base/memory.cpp +++ b/Base/memory.cpp @@ -67,90 +67,8 @@ static Memory_Table* memtab = NULL; static long memory_list_size = 0; static long memory_alloced = 0; static long memory_max_alloced = 0; -extern long MemoryTestInTimeSum(void); #endif -/* Definitionen */ - -/************************************************************************** - ROCKFLOW - Funktion: InitMemoryTest - - Aufgabe: - Oeffnet Protokolldatei "memtest.log". - - Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - - void - - - Ergebnis: - 0 bei Fehler, sonst 1 - - Programmaenderungen: - 12/1994 MSR Erste Version - 5/1997 C.Thorenz Erste Version des "in-time" Speichertests - -**************************************************************************/ -int InitMemoryTest(void) -{ -#ifdef MEMORY_TEST_IN_TIME - /* Sicherstellen, dass auch beim ersten Aufruf schon Speicher zur - Verfuegung steht */ - memtab = (Memory_Table*)malloc(sizeof(Memory_Table)); - memory_list_size = 0; - memory_alloced = 0; - memory_max_alloced = 0; -#endif - - return 1; -} - -/************************************************************************** - ROCKFLOW - Funktion: StopMemoryTest - - Aufgabe: - Schliesst Protokolldatei "memtest.log". - - Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - - void - - - Ergebnis: - - void - - - Programmaenderungen: - 12/1994 MSR Erste Version - 5/1997 C.Thorenz Erste Version des "in-time" Speichertests - -**************************************************************************/ -void StopMemoryTest(void) -{ -#ifdef MEMORY_TEST_IN_TIME -#ifdef MEMORY_STR - long i; -#endif - - if (memtab) - { - printf("\n"); -#ifdef MEMORY_STR - /* Fehler in der Speicherfreigabe anzeigen */ - i = 0; - while (i < memory_list_size) - { - if (memtab[i].address > 0) - printf("Nicht freigegeben wurde: %s Zeile %d\n", memtab[i].file, memtab[i].line); - i = i + 1; - } -#endif - /* Speicher wieder freigeben */ - free(memtab); - memtab = NULL; - - printf("\nSpeicherbilanz:\n"); - printf("Nicht freigegeben wurden %ld Bytes.\n", memory_alloced); - printf("Maximal allokiert waren %ld Bytes.\n\n", memory_max_alloced); - } -#endif -} - /************************************************************************** ROCKFLOW - Funktion: Malloc (MAlloc) @@ -464,103 +382,3 @@ void* REalloc(void* block, long bytes, char* datei, int zeile) return address; } - -#ifdef MEMORY_TEST_IN_TIME -/**************************************************************************/ -/* ROCKFLOW - Funktion: MEMORY_TEST_IN_TIME_SUM - */ -/* Aufgabe: - Ermittelt bei gesetztem MEMORY_TEST_IN_TIME - die Summe des per Malloc/Free/Realloc belegten Speichers. - */ -/* Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E void : */ -/* Ergebnis: - long summe: Belegter Speicher auf dem Heap - */ -/* Programmaenderungen: - 5/1997 C.Thorenz Erste Version des "in-time" Speichertest - */ -/**************************************************************************/ -long MemoryTestInTimeSum(void) -{ - long i, summe; - - i = 0; - summe = 0; - if (memtab) - - while (i < memory_list_size) - { - summe = summe + memtab[i].size; - i = i + 1; - } - return summe; -} -#endif -#ifdef MEMORY_TEST_IN_TIME -/**************************************************************************/ -/* ROCKFLOW - Funktion: MEMORY_TEST_IN_TIME_MEMORY - */ -/* Aufgabe: - Ermittelt bei gesetztem MEMORY_TEST_IN_TIME - die Groesse des durch den Memorytest belegten Speichers. - */ -/* Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E void : */ -/* Ergebnis: - long summe: Belegter Speicher auf dem Heap - */ -/* Programmaenderungen: - 5/1997 C.Thorenz Erste Version des "in-time" Speichertest - */ -/**************************************************************************/ -long MemoryTestInTimeMemory(void) -{ - long summe; - - summe = memory_list_size * sizeof(Memory_Table); - - return summe; -} -#endif - -#ifdef MEMORY_TEST_IN_TIME -/**************************************************************************/ -/* ROCKFLOW - Funktion: MEMORY_TEST_IN_TIME_AREA - */ -/* Aufgabe: - Ermittelt bei gesetztem MEMORY_TEST_IN_TIME - die Groesse des durch (Groesste Adresse) - (Kleinste Adresse) - angegebenen Bereichs */ -/* Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E void : */ -/* Ergebnis: - long summe: Belegter Speicher auf dem Heap - */ -/* Programmaenderungen: - 5/1997 C.Thorenz Erste Version des "in-time" Speichertest - */ -/**************************************************************************/ -long MemoryTestInTimeArea(void) -{ - long i, min_address = 0, max_address = 0; - - i = 0; - if (memtab) - { - min_address = memtab[i].address; - max_address = memtab[i].address + memtab[i].size; - - while (i < memory_list_size) - { - if (min_address > memtab[i].address) - min_address = memtab[i].address; - if (max_address < (memtab[i].address + memtab[i].size)) - max_address = memtab[i].address + memtab[i].size; - i = i + 1; - } - } - return max_address - min_address; -} -#endif diff --git a/Base/memory.h b/Base/memory.h index 20f1911b4..b69eb42c7 100644 --- a/Base/memory.h +++ b/Base/memory.h @@ -29,19 +29,6 @@ /* Andere oeffentlich benutzte Module */ -/* Deklarationen */ - -extern int InitMemoryTest(void); -/* Startet Memory-Test, 0 bei Fehler */ -extern void StopMemoryTest(void); -/* Beendet Memory-Test */ -extern long MemoryTestInTimeSum(void); -/* Gibt Summe des aktuell allokierten Speichers zurueck */ -extern long MemoryTestInTimeMemory(void); -/* Gibt Summe des Speichers fuer den Eigenbedarf zurueck */ -extern long MemoryTestInTimeArea(void); -/* Gibt Groesse des aktuell durch Addressen eingegrenzten Bereichs zurueck */ - #ifndef MEMORY_STR extern void* Malloc(long bytes); /* Ersetzt malloc */ diff --git a/OGS/rf.cpp b/OGS/rf.cpp index 785ccf720..2c21680cb 100644 --- a/OGS/rf.cpp +++ b/OGS/rf.cpp @@ -264,34 +264,6 @@ int main(int argc, char* argv[]) #endif DisplayStartMsg(); - /* Speicherverwaltung initialisieren */ - if (!InitMemoryTest()) - { - DisplayErrorMsg("Fehler: Speicherprotokoll kann nicht erstellt werden!"); - DisplayErrorMsg(" Programm vorzeitig beendet!"); - return 1; // LB changed from 0 to 1 because 0 is indicating success - } - if (argc == 1) // interactive mode - - dateiname = ReadString(); - else // non-interactive mode - { - if (argc == 2) // a model root was supplied - { - dateiname = (char*)Malloc((int)strlen(argv[1]) + 1); - dateiname = strcpy(dateiname, argv[1]); - } - else // several args supplied - if (modelRoot != "") - { - dateiname = (char*)Malloc((int)modelRoot.size() + 1); - dateiname = strcpy(dateiname, modelRoot.c_str()); - } - DisplayMsgLn(dateiname); - } - // WW DisplayMsgLn(""); - // WW DisplayMsgLn(""); - // ----------23.02.2009. WW----------------- // LB Check if file exists std::string tmpFilename = dateiname; From bbb93f50b46b756b94736577f7d45ea1dd5d4422 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Tue, 29 May 2018 13:13:22 +0200 Subject: [PATCH 6/9] Removed unused source code in some other files --- FEM/Output.cpp | 11 - FEM/files0.cpp | 29 -- FEM/files0.h | 1 - FEM/problem.cpp | 22 +- FEM/problem.h | 2 - FEM/rf_mfp_new.cpp | 405 ------------------------- FEM/rf_mmp_new.cpp | 150 ---------- FEM/rf_pcs.cpp | 726 --------------------------------------------- FEM/rf_pcs.h | 56 +--- FEM/rfmat_cp.cpp | 33 --- 10 files changed, 5 insertions(+), 1430 deletions(-) diff --git a/FEM/Output.cpp b/FEM/Output.cpp index 397be5e66..62fffdaab 100644 --- a/FEM/Output.cpp +++ b/FEM/Output.cpp @@ -1460,17 +1460,6 @@ void COutput::WriteELEValuesTECData(fstream& tec_file) tec_file << m_pcs_2->GetElementValue(i, ele_value_index_vector[j]) << " "; } } - /* - int j; - int eidx; - char ele_value_char[20]; - int no_ele_values = (int)ele_value_vector.size(); - for(j=0;j 0) + // The body of GetRFProcessProcessingAndActivation just returns 0. + // Therefore, the following destruction never happens under + // if (GetRFProcessProcessingAndActivation("MT") && GetRFProcessNumComponents() > 0) + if (GetRFProcessNumComponents() > 0) { DestroyREACT(); // SB cp_vec.clear(); // Destroy component properties vector @@ -4436,22 +4438,6 @@ void Problem::createShapeFunctionPool() } } -/************************************************************************** - FEMLib-Method: - 06/2009 OK Implementation -**************************************************************************/ -bool MODCreate() -{ - PCSConfig(); // OK - if (!PCSCheck()) // OK - { - std::cout << "Not enough data for MOD creation.\n"; - return false; - } - else - return true; -} - #ifdef BRNS // BRNS-Coupling: For writing spatially resolved reaction rates at the final iteration, diff --git a/FEM/problem.h b/FEM/problem.h index 165cdf206..c1e04356d 100644 --- a/FEM/problem.h +++ b/FEM/problem.h @@ -181,6 +181,4 @@ class Problem static const size_t max_processes = 16; }; - -extern bool MODCreate(); // OK #endif diff --git a/FEM/rf_mfp_new.cpp b/FEM/rf_mfp_new.cpp index a3a71fbdd..244a166b8 100644 --- a/FEM/rf_mfp_new.cpp +++ b/FEM/rf_mfp_new.cpp @@ -2184,75 +2184,6 @@ double MFPCalcFluidsHeatCapacity(CFiniteElementStd* assem) return heat_capacity_fluids; } -/************************************************************************** - FEMLib-Method: - Task: Master calc function - Programing: - 08/2004 OK MFP implementation based on MATCalcFluidHeatCapacity (OK) - 10/2005 YD/OK: general concept for heat capacity - overloaded function, see above, taken out, CMCD -**************************************************************************/ -/*double MFPCalcFluidsHeatCapacity(double temperature, CFiniteElementStd* assem) - { - double saturation_phase; - double heat_capacity_fluids=0.0; - int nidx0,nidx1; - //-------------------------------------------------------------------- - // MMP medium properties - bool New = false; // To be removed. WW - if(fem_msh_vector.size()>0) New = true; - //---------------------------------------------------------------------- - CFluidProperties *m_mfp = NULL; - int no_fluids =(int)mfp_vector.size(); - switch(no_fluids){ - //.................................................................... - case 1: - //m_mfp = mfp_vector[0]; - if(New) //WW - { - heat_capacity_fluids = assem->FluidProp->Density() \ - * assem->FluidProp->SpecificHeatCapacity(); - } - else - heat_capacity_fluids = assem->FluidProp->Density() \ - * assem->FluidProp->SpecificHeatCapacity(); - - break; - //.................................................................... - case 2: - if(New) //WW - { - nidx0 = GetNodeValueIndex("SATURATION1"); - nidx1 = nidx0+1; - saturation_phase = (1.-assem->pcs->m_num->ls_theta)*assem->interpolate(nidx0) - + assem->pcs->m_num->ls_theta*assem->interpolate(nidx1); - } - else - { - nidx0 = PCSGetNODValueIndex("SATURATION1",0); - nidx1 = PCSGetNODValueIndex("SATURATION1",1); - saturation_phase = (1.-assem->pcs->m_num->ls_theta)*assem->interpolate(nidx0) \ - + assem->pcs->m_num->ls_theta*assem->interpolate(nidx1); - } - m_mfp = mfp_vector[0]; - heat_capacity_fluids = saturation_phase \ - * assem->FluidProp->Density() \ - * assem->FluidProp->SpecificHeatCapacity(); - m_mfp = mfp_vector[1]; - heat_capacity_fluids += (1.0-saturation_phase) \ - * assem->FluidProp->Density() \ - * assem->FluidProp->SpecificHeatCapacity(); - break; - //.................................................................... - case 3: // Entropy based - break; - //.................................................................... - default: - cout << "Error in MFPCalcFluidsHeatCapacity: no fluid phase data" << "\n"; - } - return heat_capacity_fluids; - }*/ - /************************************************************************** FEMLib-Method: Task: Master calc function @@ -2415,342 +2346,6 @@ double MFPCalcFluidsHeatConductivity(long index, double* gp, double theta, CFini return heat_conductivity_fluids; } -//////////////////////////////////////////////////////////////////////////// -// Fluid phase change properties -#ifdef obsolete // WW -/************************************************************************** - FEMLib-Method: - Task: Vapour pressure from table - Programing: - 03/2002 OK/JdJ Implementation - last modification: -**************************************************************************/ -double MFPCalcVapourPressure(double temperature) -{ - double temperature_F; - double pressure; - double quality; - double weight; - double enthalpy; - double entropy; - double saturation_temperature; - double saturation_pressure; - double degrees_superheat; - double degrees_subcooling; - double viscosity; - double critical_velocity; - int action = 0; - double pressure_vapour; - - /* - double vapour_pressure; - double vapour_enthalpy = 2258.0; kJ/kg - double potenz; - potenz = ((1./temperature_ref)-(1./(*temperature))) * \ - ((vapour_enthalpy*FluidConstant::ComponentMolarMassWater())/FluidConstant::GasConstant()); - vapour_pressure = pressure_ref * exp(potenz); - */ - pressure = 1.e-3; /*Vorgabe eines vernueftigen Wertes*/ - pressure *= PA2PSI; /* Umrechnung Pa in psia */ - temperature -= 273.15; /* Kelvin -> Celsius */ - temperature_F = temperature * 1.8 + 32.; /* Umrechnung Celsius in Fahrenheit */ - steam67(&temperature_F, - &pressure, - &quality, - &weight, - &enthalpy, - &entropy, - &saturation_temperature, - &saturation_pressure, - °rees_superheat, - °rees_subcooling, - &viscosity, - &critical_velocity, - action); - pressure_vapour = saturation_pressure * PSI2PA; /* Umrechnung psia in Pa */ - /* - density_vapour = 0.062427962/weight; Dichte in kg/m^3 - enthalpy_vapour = enthalpy * 1055. / 0.454; Umrechnung von btu/lbm in J/kg - */ - return pressure_vapour; -} - -/************************************************************************** - FEMLib-Method: - Task: Get phase and species related enthalpies - Programing: - 03/2002 OK/JdJ Implementation - 08/2004 OK MFP implementation - last modification: -**************************************************************************/ -double CFluidProperties::Enthalpy(int comp, double temperature) -{ - double temperature_F; - double pressure; - double quality; - double weight; - double enthalpy = 0.0; - double entropy; - double saturation_temperature; - double saturation_pressure; - double degrees_superheat; - double degrees_subcooling; - double viscosity; - double critical_velocity; - int action = 0; - - if ((phase == 0) && (comp == 0)) - enthalpy = 733.0 * temperature - + (FluidConstant::GasConstant() * (temperature + 0.0)) / FluidConstant::ComponentMolarMassAir(); - else if ((phase == 0) && (comp == 1)) /* h_w^g: water species in gaseous phase */ - { - pressure = 1.e-3; /*Vorgabe eines vernuenftigen Wertes */ - pressure *= PA2PSI; /* Umrechnung Pa in psia */ - temperature -= 273.15; /* Kelvin -> Celsius */ - temperature_F = temperature * 1.8 + 32.; /* Umrechnung Celsius in Fahrenheit*/ - steam67(&temperature_F, - &pressure, - &quality, - &weight, - &enthalpy, - &entropy, - &saturation_temperature, - &saturation_pressure, - °rees_superheat, - °rees_subcooling, - &viscosity, - &critical_velocity, - action); - enthalpy = enthalpy * 1055. / 0.454; /* Umrechnung von btu/lbm in J/kg */ - } - else if ((phase == 1) && (comp == 0)) /* h_a^l: air species in liquid phase */ - { - } - else if ((phase == 1) && (comp == 1)) /* h_w^l: water species in liquid phase */ - { - } - return enthalpy; -} - -/************************************************************************** - FEMLib-Method: - Task: Get phase and species related enthalpies - Programing: - 03/2002 OK/JdJ Implementation - 08/2004 OK MFP implementation - last modification: -**************************************************************************/ -double CFluidProperties::EnthalpyPhase(long number, int comp, double* gp, double theta) -{ - double temperature; - double enthalpy = 0.0; - double mass_fraction_air, enthalpy_air, mass_fraction_water, enthalpy_water; - int nidx0, nidx1; - - nidx0 = PCSGetNODValueIndex("TEMPERATURE1", 0); - nidx1 = PCSGetNODValueIndex("TEMPERATURE1", 1); - temperature = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - - if (phase == 0) - { - if (comp < 0) - { - comp = 0; - mass_fraction_air = MassFraction(number, comp, gp, theta); - comp = 1; - mass_fraction_water = MassFraction(number, comp, gp, theta); - comp = 0; - enthalpy_air = Enthalpy(comp, temperature); - comp = 1; - enthalpy_water = Enthalpy(comp, temperature); - - enthalpy = mass_fraction_air * enthalpy_air + mass_fraction_water * enthalpy_water; - } - else if (comp >= 0) - enthalpy = Enthalpy(comp, temperature); - } - else if (phase == 1) - // - enthalpy = SpecificHeatCapacity() * temperature; - - return enthalpy; -} - -/************************************************************************** - FEMLib-Method: - Task: - Calculate Henry constant - Programing: - 02/2002 JdJ First implementation - 04/2003 JdJ temperature conversion from celsius to kelvin - last modification: -**************************************************************************/ -double MFPCalcHenryConstant(double temperature) -{ - double HenryConstant; - HenryConstant = (0.8942 + 1.47 * exp(-0.04394 * (temperature - 273.14))) * 0.0000000001; - return HenryConstant; -} - -/************************************************************************** - FEMLib-Method: - Task: - Calculate mass fractions - of gas phase X^g_a, X^g_w according to Claudius-Clapeyron law - of liquid phase X^l_a, X^l_w according to Henry law - Programing: - 01/2002 OK/JdJ First implementation - 03/2002 JdJ Gas phase in argument for pressure in mass fraction. - 04/2003 JdJ Dichteberechnung (setzt voraus, das Phase 0 gas ist). - 04/2003 JdJ Druckberechnung (setzt voraus, das Phase 0 gas ist). - 04/2004 JdJ Bugfix Gauss Points - 08/2004 OK MFP implementation - last modification: -**************************************************************************/ -double CFluidProperties::MassFraction(long number, int comp, double* gp, double theta, CFiniteElementStd* assem) -{ - double mass_fraction = 0.0; - double mass_fraction_air_in_gas, mass_fraction_air_in_liquid; - double gas_density = 0.0; - double vapour_pressure; - double temperature; - double henry_constant; - int nidx0, nidx1; - double gas_pressure; - /*--------------------------------------------------------------------------*/ - /* Get and calc independent variables */ - nidx0 = PCSGetNODValueIndex("PRESSURE1", 0); - nidx1 = PCSGetNODValueIndex("PRESSURE1", 1); - if (mode == 0) // Gauss point values - - gas_pressure = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - else // Node values - - gas_pressure = (1. - theta) * GetNodeVal(number, nidx0) + theta * GetNodeVal(number, nidx1); - nidx0 = PCSGetNODValueIndex("TEMPERATURE1", 0); - nidx1 = PCSGetNODValueIndex("TEMPERATURE1", 1); - if (mode == 0) // Gauss point values - - temperature = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - else // Node values - - temperature = (1. - theta) * GetNodeVal(number, nidx0) + theta * GetNodeVal(number, nidx1); - gas_density = Density(); - vapour_pressure = MFPCalcVapourPressure(temperature); - /*--------------------------------------------------------------------------*/ - /* Calc mass fractions */ - switch (phase) - { - case 0: /* gas phase */ - mass_fraction_air_in_gas = ((gas_pressure - vapour_pressure) * FluidConstant::ComponentMolarMassAir()) - / (FluidConstant::GasConstant() * (temperature + 0.0) * gas_density); - mass_fraction_air_in_gas = MRange(0.0, mass_fraction_air_in_gas, 1.0); - if (comp == 0) /* air specie */ - mass_fraction = mass_fraction_air_in_gas; - if (comp == 1) /* water specie */ - mass_fraction = 1.0 - mass_fraction_air_in_gas; - break; - case 1: /* liquid phase */ - henry_constant = MFPCalcHenryConstant(temperature); - mass_fraction_air_in_liquid = FluidConstant::ComponentMolarMassAir() - / (FluidConstant::ComponentMolarMassAir() - - FluidConstant::ComponentMolarMassWater() - * (1.0 - 1.0 / (henry_constant * (gas_pressure - vapour_pressure)))); - mass_fraction_air_in_liquid = MRange(0.0, mass_fraction_air_in_liquid, 1.0); - if (comp == 0) /* air specie */ - mass_fraction = mass_fraction_air_in_liquid; - if (comp == 1) /* water specie X_w^l = 1-X_a^l */ - mass_fraction = 1.0 - mass_fraction_air_in_liquid; - break; - } - /*--------------------------------------------------------------------------*/ - return mass_fraction; -} - -/************************************************************************** - FEMLib-Method: - Task: - Calculate Henry constant - Programing: - 02/2002 OK/JdJ Implementation - 08/2004 OK MFP implementation - last modification: -**************************************************************************/ -double CFluidProperties::InternalEnergy(long number, double* gp, double theta) -{ - double energy = 0.0; - int nidx0, nidx1; - double temperature, pressure; - - nidx0 = PCSGetNODValueIndex("TEMPERATURE1", 0); - nidx1 = PCSGetNODValueIndex("TEMPERATURE1", 1); - temperature = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - energy = SpecificHeatCapacity() * temperature; // YD - // energy = HeatCapacity(temperature,m_ele_fem_std) * temperature; - - // if(name.find("GAS")!=string::npos){ - if (phase == 0) - { - nidx0 = PCSGetNODValueIndex("PRESSURE1", 0); - nidx1 = PCSGetNODValueIndex("PRESSURE1", 1); - pressure = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - energy += pressure / Density(); // YD - } - return energy; -} - -/************************************************************************** - FEMLib-Method: - Task: - Calculate Henry constant - temperature in Kelvin - Programing: - 01/2003 OK/JdJ Implementation - 08/2004 OK MFP implementation - last modification: -**************************************************************************/ -double CFluidProperties::DensityTemperatureDependence(long number, int comp, double* gp, double theta) -{ - double vapour_pressure; - double dvapour_pressure_dT; - double drho_dT; - double temperature; - int nidx0, nidx1; - //---------------------------------------------------------------------- - // State functions - nidx0 = PCSGetNODValueIndex("TEMPERATURE1", 0); - nidx1 = PCSGetNODValueIndex("TEMPERATURE1", 1); - if (mode == 0) // Gauss point values - - temperature = (1. - theta) * InterpolValue(number, nidx0, gp[0], gp[1], gp[2]) - + theta * InterpolValue(number, nidx1, gp[0], gp[1], gp[2]); - else // Node values - - temperature = (1. - theta) * GetNodeVal(number, nidx0) + theta * GetNodeVal(number, nidx1); - //---------------------------------------------------------------------- - // Vapour - vapour_pressure = MFPCalcVapourPressure(temperature); - dvapour_pressure_dT = FluidConstant::ComponentMolarMassWater() * Enthalpy(comp, temperature) - / (FluidConstant::GasConstant() * temperature * temperature) * vapour_pressure; - - drho_dT = -1.0 * FluidConstant::ComponentMolarMassWater() / (FluidConstant::GasConstant() * temperature) - * (dvapour_pressure_dT - vapour_pressure / temperature); - //---------------------------------------------------------------------- - // Test - if ((phase > 0) || (comp == 0)) - { - DisplayMsgLn("MATCalcFluidDensityTemperatureDependence: Incorrect use of function"); - abort(); - } - return drho_dT; -} -#endif // if define obsolete. WW - double CFluidProperties::LiquidViscosity_expo(double T) { return viscosity0 * exp(-(T - T_0) / viscosity_T_star); diff --git a/FEM/rf_mmp_new.cpp b/FEM/rf_mmp_new.cpp index 9b5928bc6..5d29cfa3a 100644 --- a/FEM/rf_mmp_new.cpp +++ b/FEM/rf_mmp_new.cpp @@ -4148,7 +4148,6 @@ double CMediumProperties::PorosityEffectiveConstrainedSwellingConstantIonicStren /* porosity_IL = porosity_IL*satu; */ porosity_IL = porosity_max * pow(satu, beta); - // ElSetElementVal(index,PCSGetELEValueIndex("POROSITY_IL"),porosity_IL); m_pcs->SetElementValue(index, m_pcs->GetElementValueIndex("POROSITY_IL") + 1, porosity_IL); /*-----------Effective porosity calculation------------------*/ @@ -4158,7 +4157,6 @@ double CMediumProperties::PorosityEffectiveConstrainedSwellingConstantIonicStren if (porosity < porosity_min) porosity = porosity_min; - // ElSetElementVal(index,PCSGetELEValueIndex("POROSITY"),porosity); m_pcs->SetElementValue(index, m_pcs->GetElementValueIndex("POROSITY") + 1, porosity); /*-----------Swelling strain rate, xie, wang and Kolditz, (23)------------------*/ @@ -6622,78 +6620,6 @@ double CMediumProperties::PorosityEffectiveConstrainedSwelling(long index, return porosity; } -/**************************************************************************/ -/* ROCKFLOW - Funktion: CalculateSoilPorosityMethod1 - */ -/* Aufgabe: - Berechnet die Porositaet in Abhaengigkeit von der Konzentration - (Salzloesungsmodell) - */ -/* Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E SOIL_PROPERTIES *sp: Zeiger auf eine Instanz vom Typ SOIL_PROPERTIES. - */ -/* Ergebnis: - - s.o. - - */ -/* Programmaenderungen: - 02/2000 RK Erste Version - 11/2001 AH Warnung entfernt - - */ -/**************************************************************************/ -/* - double CalculateSoilPorosityMethod1(SOIL_PROPERTIES * sp, long index) - { - // double grainbound_limit = 1.73e-8; - double porosity_limit = 0.01; - int porosity_dependence_model; - double val0, val1; - static double porosity = 0.0; - static double rho; - static double theta; - static double porosity_n, density_rock; - static double dissolution_rate, solubility_coefficient; - int timelevel = 1; - - porosity_dependence_model = get_sp_porosity_dependence_model(sp); - switch (porosity_dependence_model) { - case 1: //Salzloesungsmodell - val0 = get_sp_porosity_model_field_value(sp, 0); - val1 = get_sp_porosity_model_field_value(sp, 1); - - theta = GetNumericalTimeCollocation("TRANSPORT"); - density_rock = GetSolidDensity(index); - dissolution_rate = GetTracerDissolutionRate(index, 0, 0); - porosity_n = GetElementSoilPorosity(index); - rho = GetFluidDensity(0, index, 0., 0., 0., theta); - if (GetTracerSolubilityModel(index, 0, 0) == 1) { - //if (GetRFProcessHeatReactModel()) - solubility_coefficient = - CalcTracerSolubilityCoefficient(index,0,0,1,PCSGetNODValueIndex("CONCENTRATION1",timelevel),1,PCSGetNODValueIndex("CONCENTRATION1",timelevel)); - //else solubility_coefficient = - CalcTracerSolubilityCoefficient(index,0,0,1,PCSGetNODValueIndex("CONCENTRATION1",timelevel),0,PCSGetNODValueIndex("CONCENTRATION1",timelevel)); - } else { - solubility_coefficient = GetTracerSolubilityCoefficient(index, 0, 0); - } - - porosity = porosity_n + 2 * porosity_n * dissolution_rate * val1 * rho - * (solubility_coefficient - val0) * dt / density_rock; - - if (porosity > porosity_limit) - porosity = porosity_limit; - - break; - - default: - DisplayMsgLn("Unknown porosity dependence model!"); - break; - - } - - return porosity; - } - */ - /************************************************************************** FEMLib-Method: Task: @@ -6885,34 +6811,6 @@ double CMediumProperties::TortuosityFunction(long number, double* gp, double the number = number; assem = assem; - // static int nidx0,nidx1; - //---------------------------------------------------------------------- - /*OK411 - int count_nodes; - long* element_nodes = NULL; - double primary_variable[10]; //OK To Do - int i; - int no_pcs_names =(int)pcs_name_vector.size(); - for(i=0;itimelevel = 1; // always at new time level - pcs_nval_data[pcs_nval].timelevel = pcs_secondary_function_timelevel[i]; - if (pcs_nval_data[pcs_nval].timelevel == 1) - pcs_nval_data[pcs_nval].speichern = 1; - else - pcs_nval_data[pcs_nval].speichern = 0; - pcs_nval_data[pcs_nval].laden = 0; - pcs_nval_data[pcs_nval].restart = 1; - pcs_nval_data[pcs_nval].adapt_interpol = 1; - pcs_nval_data[pcs_nval].vorgabe = 0.0; -#ifdef PCS_NOD - pcs_nval_data[pcs_nval].nval_index = pcs_nval; -#else - pcs_nval_data[pcs_nval].nval_index = anz_nval + pcs_nval; -#endif - if (type == 4 || type == 41) - if (dm_number_of_primary_nvals == 2 || (dm_number_of_primary_nvals == 3 && this->type == 41)) - { - // Block: - // STRESS_ZX, STRESS_YZ, STRAIN_ZX, STRAIN_ZY and LUMPED_STRESS - if (i == 3 || i == 4 || i == 9 || i == 10 || i == 13) - pcs_nval_data[pcs_nval].speichern = 0; - if (!problem_2d_plane_dm) - // Block STRESS_ZZ and STRAIN_ZZ - if (i == 5 || i == 11) - pcs_nval_data[pcs_nval].speichern = 0; - // if(!pcs_plasticity) - // { - // if(i==12) pcs_nval_data[pcs_nval].speichern = 0; //STRAIN_PLS - // } - } - pcs_nval++; - } -} - -/************************************************************************* - ROCKFLOW - Function: CRFProcess::PCSConfigNODValues - Task: Config node values - Programming: 02/2003 OK Implementation - last modified: - **************************************************************************/ -void CRFProcess::ConfigNODValues2(void) -{ - int i; - - number_of_nvals = 2 * GetPrimaryVNumber() + pcs_number_of_secondary_nvals; - for (i = 0; i < number_of_nvals; i++) - ModelsAddNodeValInfoStructure(pcs_nval_data[i].name, pcs_nval_data[i].einheit, pcs_nval_data[i].speichern, - pcs_nval_data[i].laden, pcs_nval_data[i].restart, pcs_nval_data[i].adapt_interpol, - pcs_nval_data[i].vorgabe); -} -#endif //#ifndef NEW_EQS //WW. 07.11.2008 -/************************************************************************** - FEMLib-Method: - Task: - Programing: - 07/2004 OK Implementation -**************************************************************************/ -int PCSGetNODValueIndex(const string& name, int timelevel) -{ - // PCS primary variables - int pcs_vector_size = (int)pcs_vector.size(); - int i, j; - CRFProcess* m_pcs = NULL; - if (pcs_vector_size > 0) - for (i = 0; i < pcs_vector_size; i++) - { - m_pcs = pcs_vector[i]; - for (j = 0; j < m_pcs->number_of_nvals; j++) - if ((name.compare(m_pcs->pcs_nval_data[j].name) == 0) - && (m_pcs->pcs_nval_data[j].timelevel == timelevel)) - return m_pcs->pcs_nval_data[j].nval_index; - } - cout << "Error in PCSGetNODValueIndex: " << name << "\n"; - return -1; -} - -////////////////////////////////////////////////////////////////////////// -// Configuration ELE -////////////////////////////////////////////////////////////////////////// - -/************************************************************************* - ROCKFLOW - Function: - Task: Config element values - Programming: 02/2003 OK Implementation - last modified: - 06/2004 WW - **************************************************************************/ -void CRFProcess::ConfigELEValues1(void) -{ - int i; - if (pcs_number_of_evals) - pcs_eval_data = (PCS_EVAL_DATA*)Malloc(pcs_number_of_evals * sizeof(PCS_EVAL_DATA)); - for (i = 0; i < pcs_number_of_evals; i++) - { - // pcs_eval_data[i] = (PCS_EVAL_DATA *) Malloc(sizeof(PCS_EVAL_DATA)); - strcpy(pcs_eval_data[i].name, pcs_eval_name[i]); - strcpy(pcs_eval_data[i].einheit, pcs_eval_unit[i]); - pcs_eval_data[i].speichern = 1; - pcs_eval_data[i].laden = 0; - pcs_eval_data[i].restart = 1; - pcs_eval_data[i].adapt_interpol = 1; - pcs_eval_data[i].vorgabe = 0.0; - pcs_eval_data[i].index = anz_eval + i; - pcs_eval_data[i].eval_index = anz_eval + i; // SB - } -} - -/************************************************************************* - ROCKFLOW - Function: - Task: Config element values - Programming: 04/2003 OK Implementation - last modified: - **************************************************************************/ -void CRFProcess::ConfigELEValues2(void) -{ - int i; - for (i = 0; i < pcs_number_of_evals; i++) - ModelsAddElementValInfoStructure(pcs_eval_data[i].name, pcs_eval_data[i].einheit, pcs_eval_data[i].speichern, - pcs_eval_data[i].laden, pcs_eval_data[i].restart, - pcs_eval_data[i].adapt_interpol, pcs_eval_data[i].vorgabe); -} - -/************************************************************************* - ROCKFLOW - Function: CRFProcess::PCSGetELEValueIndex - Task: Provide index for element values - Programming: 08/2003 SB Implementation - last modified: - **************************************************************************/ -int PCSGetELEValueIndex(char* name) -{ - int i; - CRFProcess* m_process = NULL; - int j; - int no_processes = (int)pcs_vector.size(); - for (j = 0; j < no_processes; j++) - { - m_process = pcs_vector[j]; - for (i = 0; i < m_process->pcs_number_of_evals; i++) - if (strcmp(m_process->pcs_eval_data[i].name, name) == 0) - return m_process->pcs_eval_data[i].eval_index; - } - printf("PCSGetELEValueIndex Alert\n"); - printf("%s \n", name); - return -1; -} - -////////////////////////////////////////////////////////////////////////// -// Configuration ELE matrices -////////////////////////////////////////////////////////////////////////// - /************************************************************************** FEMLib-Method: Task: Activate or deactivate elements specified in .pcs file @@ -8279,49 +8048,6 @@ int CRFProcess::GetNODValueIndex(const string& name, int timelevel) return -1; } -/////////////////////////////////////////////////////////////////////////// -// Specials -/////////////////////////////////////////////////////////////////////////// - -/*------------------------------------------------------------------------- - ROCKFLOW - Function: PCSRestart - Task: Insert process to list - Programming: - 06/2003 OK Implementation - 11/2004 OK file_name_base - last modified: - -------------------------------------------------------------------------*/ -void PCSRestart() -{ - /*OK411 - int j; - CRFProcess *m_pcs = NULL; - int nidx0,nidx1; - int i; - int no_processes =(int)pcs_vector.size(); - if(no_processes==0) - return; //OK41 - int ok = 0; - //---------------------------------------------------------------------- - string file_name_base = pcs_vector[0]->file_name_base; - //OK ok = ReadRFRRestartData(file_name_base); - if(ok==0){ - cout << "RFR: no restart data" << "\n"; - return; - } - //---------------------------------------------------------------------- - for(i=0;iGetPrimaryVNumber();j++) { - // timelevel=0; - nidx0 = m_pcs->GetNodeValueIndex(m_pcs->GetPrimaryVName(j)); - // timelevel= 1; - nidx1 = nidx0+1; - CopyNodeVals(nidx1,nidx0); - } - } - */ -} /************************************************************************** FEMLib-Method: @@ -8372,226 +8098,6 @@ void CRFProcess::PCSMoveNOD(void) } } -/************************************************************************** - FEMLib-Method: - Task: - Programing: - 09/2004 OK Implementation - 10/2004 OK 2nd version -**************************************************************************/ -std::string PCSProblemType() -{ - std::string pcs_problem_type; - size_t no_processes(pcs_vector.size()); - - for (size_t i = 0; i < no_processes; i++) - { - switch (pcs_vector[i]->getProcessType()) - { - case FiniteElement::LIQUID_FLOW: - pcs_problem_type = "LIQUID_FLOW"; - break; - case FiniteElement::OVERLAND_FLOW: - pcs_problem_type = "OVERLAND_FLOW"; - break; - case FiniteElement::GROUNDWATER_FLOW: - pcs_problem_type = "GROUNDWATER_FLOW"; - break; - case FiniteElement::TWO_PHASE_FLOW: - pcs_problem_type = "TWO_PHASE_FLOW"; - break; - case FiniteElement::RICHARDS_FLOW: // MX test 04.2005 - pcs_problem_type = "RICHARDS_FLOW"; - break; - case FiniteElement::DEFORMATION: - if (pcs_problem_type.empty()) - pcs_problem_type = "DEFORMATION"; - else - pcs_problem_type += "+DEFORMATION"; - break; - case FiniteElement::DEFORMATION_FLOW: - if (pcs_problem_type.empty()) - pcs_problem_type = "DEFORMATION"; - else - pcs_problem_type += "+DEFORMATION"; - break; - case FiniteElement::HEAT_TRANSPORT: - if (pcs_problem_type.empty()) - pcs_problem_type = "HEAT_TRANSPORT"; - else - pcs_problem_type += "+HEAT_TRANSPORT"; - break; - case FiniteElement::MASS_TRANSPORT: - if (pcs_problem_type.empty()) - pcs_problem_type = "MASS_TRANSPORT"; - else - pcs_problem_type += "+MASS_TRANSPORT"; - break; - case FiniteElement::FLUID_MOMENTUM: - if (pcs_problem_type.empty()) - pcs_problem_type = "FLUID_MOMENTUM"; - else - pcs_problem_type += "+FLUID_MOMENTUM"; - break; - case FiniteElement::RANDOM_WALK: - if (pcs_problem_type.empty()) - pcs_problem_type = "RANDOM_WALK"; - else - pcs_problem_type += "+RANDOM_WALK"; - break; - default: - pcs_problem_type = ""; - } - } - - // //---------------------------------------------------------------------- - // CRFProcess* m_pcs = NULL; - // // H process - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[0]) { - // case 'L': - // pcs_problem_type = "LIQUID_FLOW"; - // break; - // // case 'U': - // // pcs_problem_type = "UNCONFINED_FLOW"; - // // break; - // case 'O': - // pcs_problem_type = "OVERLAND_FLOW"; - // break; - // case 'G': - // pcs_problem_type = "GROUNDWATER_FLOW"; - // break; - // case 'T': - // pcs_problem_type = "TWO_PHASE_FLOW"; - // break; - // case 'C': - // pcs_problem_type = "COMPONENTAL_FLOW"; - // break; - // case 'R': //MX test 04.2005 - // pcs_problem_type = "RICHARDS_FLOW"; - // break; - // } - // } - // //---------------------------------------------------------------------- - // // M process - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[0]) { - // case 'D': - // if (pcs_problem_type.empty()) - // pcs_problem_type = "DEFORMATION"; - // else - // pcs_problem_type += "+DEFORMATION"; - // break; - // } - // } - // //---------------------------------------------------------------------- - // // T process - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[0]) { - // case 'H': - // if (pcs_problem_type.empty()) - // pcs_problem_type = "HEAT_TRANSPORT"; - // else - // pcs_problem_type += "+HEAT_TRANSPORT"; - // break; - // } - // } - // //---------------------------------------------------------------------- - // // CB process - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[0]) { - // case 'M': - // if (pcs_problem_type.empty()) - // pcs_problem_type = "MASS_TRANSPORT"; - // else - // pcs_problem_type += "+MASS_TRANSPORT"; - // break; - // } - // } - // //---------------------------------------------------------------------- - // //---------------------------------------------------------------------- - // // FM process - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[0]) { - // case 'F': - // if (pcs_problem_type.empty()) - // pcs_problem_type = "FLUID_MOMENTUM"; - // else - // pcs_problem_type += "+FLUID_MOMENTUM"; - // break; - // } - // } - // //---------------------------------------------------------------------- - // for (i = 0; i < no_processes; i++) { - // m_pcs = pcs_vector[i]; - // switch (m_pcs->pcs_type_name[7]) { // _pcs_type_name[7] should be 'W' because 'R' is reserved for Richard - //Flow. - // case 'W': - // if (pcs_problem_type.empty()) - // pcs_problem_type = "RANDOM_WALK"; - // else - // pcs_problem_type += "+RANDOM_WALK"; - // break; - // } - // } - - return pcs_problem_type; -} - -/************************************************************************** - FEMLib-Method: - Task: - Programing: - 11/2004 OK Implementation -**************************************************************************/ -// void CRFProcess::CalcELEMassFluxes(void) -//{ -/*OK411 - int i; - double e_value = -1.0; - string e_value_name; - double geo_factor, density; - double velocity = 0.0; - int e_idx; - int phase = 0; - long e; - CMediumProperties* m_mmp = NULL; - CFluidProperties* m_mfp = NULL; - m_mfp = mfp_vector[phase]; //OK ToDo - //====================================================================== - for(e=0;egeo_area; - density = m_mfp->Density(); - for(i=0;i_pcs_type_name.find("DEFORMATION") != string::npos) - if (isDeformationProcess(m_pcs->getProcessType())) - pcs_deform = true; - // if (m_pcs->_pcs_type_name.find("FLOW") != string::npos) - if (isFlowProcess(m_pcs->getProcessType())) - pcs_flow = true; - } - - if (strcmp(rfpp_type, "SD") == 0) - { - if (pcs_flow && pcs_deform) - return 1; - } - else - cout << "GetRFProcessProcessing - to be removed" - << "\n"; - return 0; -} - -int GetRFProcessProcessingAndActivation(const char*) -{ - cout << "GetRFProcessProcessingAndActivation - to be removed" - << "\n"; - return 0; -} - long GetRFProcessNumComponents(void) { // DisplayMsgLn("GetRFProcessNumComponents - to be removed"); @@ -8843,13 +8271,6 @@ long GetRFProcessNumComponents(void) return no_components; } -int GetRFControlModex(void) -{ - cout << "GetRFControlModex - to be removed" - << "\n"; - return 0; -} - int GetRFProcessDensityFlow(void) { if (show_onces_density) @@ -8859,138 +8280,6 @@ int GetRFProcessDensityFlow(void) return 0; } -int GetRFProcessNumContinua(void) -{ - cout << "GetRFProcessNumContinua - to be removed" - << "\n"; - return 0; -} - -int GetRFProcessNumElectricFields(void) -{ - cout << "GetRFProcessNumElectricFields - to be removed" - << "\n"; - return 0; -} - -int GetRFProcessNumTemperatures(void) -{ - cout << "GetRFProcessNumTemperatures - to be removed" - << "\n"; - return -1; -} - -int GetRFProcessSimulation(void) -{ - cout << "GetRFProcessSimulation - to be removed" - << "\n"; - return -1; -} - -/************************************************************************** - ROCKFLOW - Funktion: ModelsAddNodeValInfoStructure - - Aufgabe: - Fuellt die Knotendaten-Infostruktur mit den zugehoerigen Modelldaten. - Wird vom Modell der Reihe nach fuer jede Knotengroesse aufgerufen. - - Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E:char *name :Name der Knotengroesse fuer Ergebnisdatei - E:char *einheit :Name der phys. Einheit fuer Ergebnisdatei - E:int speichern :Werte sollen gespeichert werden (0/1) - E:int laden :Werte sollen geladen werden falls vorhanden (0/1) - E:int restart :Werte sollen bei Restart geladen werden (0/1) - E:int adapt_interpol :Werte sollen beim verfeinern auf Kinder interpoliert (0/1) - E:double vorgabe :Vorgabe falls keine Restartdaten oder Anfangsbedingungen vorhanden sind - - Ergebnis: - Knotenindex der gerade vergeben wurde - - Programmaenderungen: - 09/2000 CT Erste Version - -**************************************************************************/ -int ModelsAddNodeValInfoStructure(char* name, char* einheit, int speichern, int laden, int restart, int adapt_interpol, - double vorgabe) -{ - anz_nval++; - nval_data = (NvalInfo*)Realloc(nval_data, anz_nval * sizeof(NvalInfo)); - - nval_data[anz_nval - 1].name = NULL; - nval_data[anz_nval - 1].einheit = NULL; - - if (name) - { - nval_data[anz_nval - 1].name = (char*)Malloc(((int)strlen(name) + 1) * sizeof(char)); - strcpy(nval_data[anz_nval - 1].name, name); - } - if (einheit) - { - nval_data[anz_nval - 1].einheit = (char*)Malloc(((int)strlen(einheit) + 1) * sizeof(char)); - strcpy(nval_data[anz_nval - 1].einheit, einheit); - } - - nval_data[anz_nval - 1].speichern = speichern; - nval_data[anz_nval - 1].laden = laden; - nval_data[anz_nval - 1].restart = restart; - nval_data[anz_nval - 1].adapt_interpol = adapt_interpol; - nval_data[anz_nval - 1].vorgabe = vorgabe; - - return anz_nval - 1; -} - -/************************************************************************** - ROCKFLOW - Funktion: ModelsAddElementValInfoStructure - - Aufgabe: - Fuellt die Elementdaten-Infostruktur mit den zugehoerigen Modelldaten. - Wird vom Modell der Reihe nach fuer jede Elementgroesse aufgerufen. - - Formalparameter: (E: Eingabe; R: Rueckgabe; X: Beides) - E:char *name :Name der Elementgroesse fuer Ergebnisdatei - E:char *einheit :Name der phys. Einheit fuer Ergebnisdatei - E:int speichern :Werte sollen gespeichert werden (0/1) - E:int laden :Werte sollen geladen werden falls vorhanden (0/1) - E:int restart :Werte sollen bei Restart geladen werden (0/1) - E:int adapt_interpol :Werte sollen beim verfeinern auf Kinder interpoliert (0/1) - E:double vorgabe :Vorgabe falls keine Restartdaten oder Anfangsbedingungen vorhanden sind - - Ergebnis: - Elementindex der gerade vergeben wurde - - Programmaenderungen: - 09/2000 CT Erste Version - -**************************************************************************/ -int ModelsAddElementValInfoStructure(char* name, char* einheit, int speichern, int laden, int restart, - int adapt_interpol, double vorgabe) -{ - anz_eval++; - eval_data = (EvalInfo*)Realloc(eval_data, anz_eval * sizeof(EvalInfo)); - - eval_data[anz_eval - 1].name = NULL; - eval_data[anz_eval - 1].einheit = NULL; - - if (name) - { - eval_data[anz_eval - 1].name = (char*)Malloc(((int)strlen(name) + 1) * sizeof(char)); - strcpy(eval_data[anz_eval - 1].name, name); - } - if (einheit) - { - eval_data[anz_eval - 1].einheit = (char*)Malloc(((int)strlen(einheit) + 1) * sizeof(char)); - strcpy(eval_data[anz_eval - 1].einheit, einheit); - } - - eval_data[anz_eval - 1].speichern = speichern; - eval_data[anz_eval - 1].laden = laden; - eval_data[anz_eval - 1].restart = restart; - eval_data[anz_eval - 1].adapt_interpol = adapt_interpol; - eval_data[anz_eval - 1].vorgabe = vorgabe; - - return anz_eval - 1; -} - /************************************************************************** FEMLib-Method: Task: @@ -14711,21 +14000,6 @@ void CRFProcess::configMaterialParameters() } } -/************************************************************************* - GeoSys - Function: - 06/2009 OK Implementation - **************************************************************************/ -bool PCSConfig() -{ - bool some_thing_done = false; - for (int i = 0; i < (int)pcs_vector.size(); i++) // OK - { - pcs_vector[i]->Config(); - some_thing_done = true; - } - return some_thing_done; -} - /************************************************************************* GeoSys-Function: CalGPVelocitiesfromEclipse Task: Calculate gauss point velocities from Eclipse solution diff --git a/FEM/rf_pcs.h b/FEM/rf_pcs.h index 798c73771..db1557223 100644 --- a/FEM/rf_pcs.h +++ b/FEM/rf_pcs.h @@ -659,7 +659,6 @@ class CRFProcess : public ProcessInfo #if !defined(USE_PETSC) && !defined(NEW_EQS) // && defined(other parallel libs)//03~04.3012. WW void ConfigNODValues1(void); - void ConfigNODValues2(void); void CreateNODValues(void); void SetNODValues(); // OK void CalcFluxesForCoupling(); // MB @@ -679,8 +678,6 @@ class CRFProcess : public ProcessInfo void CopyCouplingNODValues(); // WWvoid CalcFluxesForCoupling(); //MB // Configuration 2 - ELE - void ConfigELEValues1(void); - void ConfigELEValues2(void); void CreateELEValues(void); void CreateELEGPValues(); void AllocateMemGPoint(); // NEW Gauss point values for CFEMSH WW @@ -984,17 +981,7 @@ extern CRFProcess* PCSGetFlow(); // OK//JT extern CRFProcess* PCSGetHeat(); // JT extern CRFProcess* PCSGetMass(size_t component_number); // JT extern CRFProcess* PCSGetDeformation(); // JT -extern bool PCSConfig(); // -// NOD -extern int PCSGetNODValueIndex(const std::string&, int); -extern double PCSGetNODValue(long, char*, int); -extern void PCSSetNODValue(long, const std::string&, double, int); -// ELE -extern int PCSGetELEValueIndex(char*); -extern double PCSGetELEValue(long index, double* gp, double theta, const std::string& nod_fct_name); -// Specials -extern void PCSRestart(); -extern std::string PCSProblemType(); + // PCS global variables extern int pcs_no_components; extern bool pcs_monolithic_flow; @@ -1017,50 +1004,9 @@ extern std::string GetPFNamebyCPName(std::string line_string); extern int memory_opt; -typedef struct /* Knotenwert-Informationen */ -{ - char* name; /* Name der Knotengroesse */ - char* einheit; /* Einheit der Knotengroesse */ - int speichern; /* s.u., wird Wert gespeichert ? */ - int laden; /* s.u., wird Wert zu Beginn geladen ? */ - int restart; /* s.u., wird Wert bei Restart geladen ? */ - int adapt_interpol; /* Soll Wert bei Adaption auf Kinder interpoliert werden? */ - double vorgabe; /* Default-Wert fuer Initialisierung */ -} NvalInfo; -extern int anz_nval; /* Anzahl der Knotenwerte */ -extern int anz_nval0; /* Anzahl der Knotenwerte */ -extern NvalInfo* nval_data; -extern int ModelsAddNodeValInfoStructure(char*, char*, int, int, int, int, double); - -typedef struct /* Elementwert-Informationen */ -{ - char* name; /* Name der Elementgroesse */ - char* einheit; /* Einheit der Elementgroesse */ - int speichern; /* s.u., wird Wert gespeichert ? */ - int laden; /* s.u., wird Wert zu Beginn geladen ? */ - int restart; /* s.u., wird Wert bei Restart geladen ? */ - int adapt_interpol; /* Soll Wert bei Adaption auf Kinder vererbt werden? */ - double vorgabe; /* Default-Wert fuer Initialisierung */ -} EvalInfo; -extern int anz_eval; /* Anzahl der Elementwerte */ -extern EvalInfo* eval_data; -extern int ModelsAddElementValInfoStructure(char*, char*, int, int, int, int, double); - -extern int GetRFControlGridAdapt(void); -extern int GetRFControlModel(void); -extern int GetRFProcessChemicalModel(void); -extern int GetRFProcessFlowModel(void); -extern int GetRFProcessHeatReactModel(void); extern int GetRFProcessNumPhases(void); -extern int GetRFProcessProcessing(char*); -extern int GetRFProcessProcessingAndActivation(const char*); extern long GetRFProcessNumComponents(void); -extern int GetRFControlModex(void); extern int GetRFProcessDensityFlow(void); -extern int GetRFProcessNumContinua(void); -extern int GetRFProcessNumElectricFields(void); -extern int GetRFProcessNumTemperatures(void); -extern int GetRFProcessSimulation(void); extern void initializeConstrainedProcesses(std::vector& pcs_vector); diff --git a/FEM/rfmat_cp.cpp b/FEM/rfmat_cp.cpp index 5fb41517c..f8dbcd5ac 100644 --- a/FEM/rfmat_cp.cpp +++ b/FEM/rfmat_cp.cpp @@ -878,39 +878,6 @@ double CompProperties::CalcDiffusionCoefficientCP(long index, double theta, CRFP case 5: case 6: case 7: -/*OK411 - { - count_nodes = ElNumberOfNodes[ElGetElementType(index) - 1]; - element_nodes = ElGetElementNodes(index); - p_ind = 1; - sprintf(name_pres,"PRESSURE%ld",transport_phase+1); - // SB:GS4 todo p_idx = m_pcs->PCSGetNODValueIndexNew(name_pres,1); - if (p_ind) { - pressure_average = 0.0; - for (i = 0; i < count_nodes; i++) - pressure_average += GetNodeVal(element_nodes[i], p_idx); - pressure_average /= count_nodes; - } - t_ind = 1; - if (! GetRFProcessHeatReactModel()){ - DisplayMsgLn("No heat transport active, no temperature available for determination of diffusion coefficient "); - break; - }; - - sprintf(name_heat,"TEMPERATURE%ld",transport_phase+1); - // SB:GS4 todo m_pcs->PCSGetNODValueIndexNew(name_heat,1); - if (t_ind) { - temperature_average = 0.0; - for (i = 0; i < count_nodes; i++) - temperature_average += GetNodeVal(element_nodes[i], t_idx); - temperature_average /= count_nodes; - } - eta = mfp_vector[transport_phase]->Viscosity(); //GetFluidViscosity(transport_phase, index, 0., 0., 0., theta); - diffusion_coefficient = CalcDiffusionCoefficientCP_Method1(index, temperature_average, pressure_average, eta); - element_nodes = NULL; - break; - } - */ #ifdef GEM_REACT case 8: /* Archies law De = Dp * tau * tau_zero* poros with tau = poros^m/poros0^m ....tortuosity = tau_zero * tau depends on porosity and tau_zero is multiplicated (like porosity) outside this function */ From 89ed00dfa718facf70804a0dcc93304f8810a9b3 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Tue, 29 May 2018 16:24:54 +0200 Subject: [PATCH 7/9] Re-written the reading of command line arguments --- FEM/files0.cpp | 2 +- FEM/problem.cpp | 4 +- FEM/problem.h | 2 +- OGS/rf.cpp | 234 +++++++++++++++++++++++++++--------------------- 4 files changed, 134 insertions(+), 108 deletions(-) diff --git a/FEM/files0.cpp b/FEM/files0.cpp index 41426c81c..e4ffc48eb 100644 --- a/FEM/files0.cpp +++ b/FEM/files0.cpp @@ -174,7 +174,7 @@ static bool checkFormatOfInputFiles(const std::string& basename) last modified: OK 16.10.2002 */ /**************************************************************************/ -int ReadData(char* dateiname, GEOLIB::GEOObjects& geo_obj, std::string& unique_name) +int ReadData(const char* dateiname, GEOLIB::GEOObjects& geo_obj, std::string& unique_name) { #if defined(USE_MPI) // WW if (myrank == 0) diff --git a/FEM/problem.cpp b/FEM/problem.cpp index 570c0f8de..938397cf5 100644 --- a/FEM/problem.cpp +++ b/FEM/problem.cpp @@ -44,7 +44,7 @@ /*------------------------------------------------------------------------*/ // Data file // OK411 -extern int ReadData(char*, GEOLIB::GEOObjects& geo_obj, std::string& unique_name); +extern int ReadData(const char*, GEOLIB::GEOObjects& geo_obj, std::string& unique_name); /* PCS */ #include "pcs_dm.h" #include "rf_pcs.h" @@ -112,7 +112,7 @@ using process::CRFProcessDeformation; PreTimeloop Modification: ***************************************************************************/ -Problem::Problem (char* filename) : +Problem::Problem (const char* filename) : dt0(0.), print_result(true), _linear_shapefunction_pool(NULL), _quadr_shapefunction_pool(NULL), _geo_obj (new GEOLIB::GEOObjects), _geo_name (filename), diff --git a/FEM/problem.h b/FEM/problem.h index c1e04356d..0b5f05953 100644 --- a/FEM/problem.h +++ b/FEM/problem.h @@ -39,7 +39,7 @@ typedef double (Problem::*ProblemMemFn)(void); class Problem { public: - Problem(char* filename = NULL); + Problem(const char* filename = NULL); ~Problem(); void Euler_TimeDiscretize(); void RosenBrock_TimeDiscretize(){}; diff --git a/OGS/rf.cpp b/OGS/rf.cpp index 2c21680cb..b7b66c643 100644 --- a/OGS/rf.cpp +++ b/OGS/rf.cpp @@ -28,6 +28,8 @@ * */ #include "BuildInfo.h" +#include + #if defined(USE_MPI) || defined(USE_MPI_PARPROC) || defined(USE_MPI_REGSOIL) || defined(USE_MPI_GEMS) \ || defined(USE_MPI_KRC) #include "par_ddc.h" @@ -104,91 +106,10 @@ double elapsed_time_mpi; /**************************************************************************/ int main(int argc, char* argv[]) { - /* parse command line arguments */ - std::string anArg; - std::string modelRoot; #if defined(USE_MPI) || defined(USE_MPI_PARPROC) || defined(USE_MPI_REGSOIL) || defined(USE_MPI_GEMS) \ || defined(USE_MPI_KRC) int nb_ddc = 0; // number of cores for DDC related processes #endif - - for (int i = 1; i < argc; i++) - { - anArg = std::string(argv[i]); - if (anArg == "--help" || anArg == "-h") - { - std::cout << "Usage: ogs [MODEL_ROOT] [OPTIONS]\n" - << "Where OPTIONS are:\n" - << " -h [--help] print this message and exit\n" - << " -b [--build-info] print build info and exit\n" - << " --output-directory DIR put output files into DIR\n" - << " --version print ogs version and exit" - << "\n"; - continue; - } - if (anArg == "--build-info" || anArg == "-b") - { - std::cout << "ogs version: " << BuildInfo::OGS_VERSION << "\n" - << "ogs date: " << BuildInfo::OGS_DATE << "\n"; - std::cout << "git commit info: " << BuildInfo::GIT_COMMIT_INFO << "\n"; - std::cout << "build timestamp: " << BuildInfo::BUILD_TIMESTAMP << "\n"; - continue; - } - if (anArg == "--version") - { - std::cout << BuildInfo::OGS_VERSION << "\n"; - continue; - } - if (anArg == "--model-root" || anArg == "-m") - { - if (i + 1 >= argc) - { - std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; - std::exit(EXIT_FAILURE); - } - modelRoot = std::string(argv[++i]); - continue; - } - if (anArg == "--output-directory") - { - if (i + 1 >= argc) - { - std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; - std::exit(EXIT_FAILURE); - } - std::string path = argv[++i]; - - if (!path.empty()) - defaultOutputPath = path; - continue; - } -#if defined(USE_MPI) || defined(USE_MPI_PARPROC) || defined(USE_MPI_REGSOIL) || defined(USE_MPI_GEMS) \ - || defined(USE_MPI_KRC) - std::string decompositions; - if (anArg == "--domain-decomposition" || anArg == "-ddc") - { - decompositions = std::string(argv[++i]); - nb_ddc = atoi(decompositions.c_str()); - continue; - } -#endif - // anything left over must be the model root, unless already found - if (modelRoot == "") - modelRoot = std::string(argv[i]); - } // end of parse argc loop - - if (argc > 1 && modelRoot == "") // non-interactive mode and no model given - exit(0); // e.g. just wanted the build info - - std::string solver_pkg_name = BuildInfo::SOLVER_PACKAGE_NAME; - // No default linear solver package is in use. - if (solver_pkg_name.find("Default") == std::string::npos) - { - std::cout << "\nWarning: " << solver_pkg_name << " other than the OGS default one is in use." << std::endl; - std::cout << " The solver setting may need to be adjusted for the solution accuracy!" << std::endl; - } - - char* dateiname(NULL); #ifdef SUPERCOMPUTER // ********************************************************************* // buffered output ... important for performance on cray @@ -249,36 +170,122 @@ int main(int argc, char* argv[]) // Initialization of the lis solver. lis_initialize(&argc, &argv); #endif -/*========================================================================*/ -/* Kommunikation mit Betriebssystem */ -/* Timer fuer Gesamtzeit starten */ -#ifdef TESTTIME - TStartTimer(0); -#endif -/* Intro ausgeben */ + + std::vector arg_strings; + if (argc > 1) + { + for (int i = 1; i < argc; i++) + { + arg_strings.push_back(std::string(argv[i])); + } + } + else + { #if defined(USE_MPI) // WW - if (myrank == 0) + if (myrank == 0) #endif #ifdef USE_PETSC if (rank == 0) #endif + DisplayStartMsg(); - DisplayStartMsg(); + std::string s_buff; + std::stringstream ss; + getline(std::cin, s_buff); + ss.str(s_buff); + while (!ss.eof()) + { + ss >> s_buff; + arg_strings.push_back(s_buff); + } + } - // LB Check if file exists - std::string tmpFilename = dateiname; - tmpFilename.append(".pcs"); - if (!IsFileExisting(tmpFilename)) + for (std::size_t i = 0; i < arg_strings.size(); i++) { - std::cout << " Error: Cannot find file " << dateiname << "\n"; - return 1; - } + const std::string anArg = arg_strings[i]; + if (anArg == "--help" || anArg == "-h") + { + std::cout << "Usage: ogs [MODEL_ROOT] [OPTIONS]\n" + << "Where OPTIONS are:\n" + << " -h [--help] print this message and exit\n" + << " -b [--build-info] print build info and exit\n" + << " --output-directory DIR put output files into DIR\n" + << " --version print ogs version and exit" + << "\n"; + exit(0); + } + if (anArg == "--build-info" || anArg == "-b") + { + std::cout << "ogs version: " << BuildInfo::OGS_VERSION << "\n" + << "ogs date: " << BuildInfo::OGS_DATE << "\n"; + std::cout << "git commit info: " << BuildInfo::GIT_COMMIT_INFO << "\n"; + std::cout << "build timestamp: " << BuildInfo::BUILD_TIMESTAMP << "\n"; + exit(0); + } + if (anArg == "--version") + { + std::cout << BuildInfo::OGS_VERSION << "\n"; + exit(0); + } + if (anArg == "--model-root" || anArg == "-m") + { + if (i + 1 >= arg_strings.size()) + { + std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; + std::exit(EXIT_FAILURE); + } +// modelRoot = std::string(arg_strings[++i]); + continue; + } + if (anArg == "--output-directory") + { + if (i + 1 >= arg_strings.size()) + { + std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; + std::exit(EXIT_FAILURE); + } + std::string path = arg_strings[++i]; - // If no option is given, output files are placed in the same directory as the input files - if (defaultOutputPath.empty()) - defaultOutputPath = pathDirname(std::string(dateiname)); + if (!path.empty()) + defaultOutputPath = path; + continue; + } +#if defined(USE_MPI) || defined(USE_MPI_PARPROC) || defined(USE_MPI_REGSOIL) || defined(USE_MPI_GEMS) \ + || defined(USE_MPI_KRC) + std::string decompositions; + if (anArg == "--domain-decomposition" || anArg == "-ddc") + { + if (i + 1 >= arg_strings.size()) + { + std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; + std::exit(EXIT_FAILURE); + } + decompositions = std::string(arg_strings[++i]); + nb_ddc = atoi(decompositions.c_str()); + continue; + } +#endif + else + { + FileName = arg_strings[i]; + } + } // end of parse argc loop + + if (argc > 1) + { +#if defined(USE_MPI) // WW + if (myrank == 0) +#endif +#ifdef USE_PETSC + if (rank == 0) +#endif + DisplayStartMsg(); + } + +#ifdef TESTTIME + TStartTimer(0); +#endif - FileName = dateiname; size_t indexChWin, indexChLinux; indexChWin = indexChLinux = 0; indexChWin = FileName.find_last_of('\\'); @@ -288,8 +295,29 @@ int main(int argc, char* argv[]) FilePath = FileName.substr(0, indexChWin) + "\\"; else if (indexChLinux != std::string::npos) FilePath = FileName.substr(0, indexChLinux) + "/"; + // If no option is given, output files are placed in the same directory as the input files + if (defaultOutputPath.empty()) + defaultOutputPath = FilePath; + + std::string solver_pkg_name = BuildInfo::SOLVER_PACKAGE_NAME; + // No default linear solver package is in use. + if (solver_pkg_name.find("Default") == std::string::npos) + { + std::cout << "\nWarning: " << solver_pkg_name << " other than the OGS default one is in use." << std::endl; + std::cout << " The solver setting may need to be adjusted for the solution accuracy!" << std::endl; + } + + // LB Check if file exists + std::string tmpFilename = FileName; + tmpFilename.append(".pcs"); + if (!IsFileExisting(tmpFilename)) + { + std::cout << " Error: Cannot find file " << FileName << "\n"; + return 1; + } + // ---------------------------WW - Problem* aproblem = new Problem(dateiname); + Problem* aproblem = new Problem(FileName.data()); #ifdef USE_PETSC aproblem->setRankandSize(rank, r_size); #endif @@ -350,8 +378,6 @@ int main(int argc, char* argv[]) #endif /*--------- LIS Finalize ------------------*/ - free(dateiname); - #ifdef USE_PETSC // kg44 quick fix to compile PETSC with version PETSCV3.4 #if (PETSC_VERSION_NUMBER > 3030) From 806b5df083572f87726ff385d7f89c3e8fd81498 Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Wed, 30 May 2018 09:49:32 +0200 Subject: [PATCH 8/9] Removed compiler related declaration in FileTools.cpp --- Base/FileTools.cpp | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/Base/FileTools.cpp b/Base/FileTools.cpp index 9add2f9b0..034c986c0 100644 --- a/Base/FileTools.cpp +++ b/Base/FileTools.cpp @@ -155,9 +155,7 @@ std::string getCwd() #ifdef WIN32 _getcwd(cwd, FILENAME_MAX); #else - char* unused __attribute__((unused)); - unused = getcwd(cwd, FILENAME_MAX); + getcwd(cwd, FILENAME_MAX); #endif - return cwd; } From 60bff175ff53584b98a7713c9e291d6d7acefd2c Mon Sep 17 00:00:00 2001 From: Wenqing Wang Date: Wed, 30 May 2018 09:50:08 +0200 Subject: [PATCH 9/9] Removed one unused command argument --- OGS/rf.cpp | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/OGS/rf.cpp b/OGS/rf.cpp index b7b66c643..bdc597352 100644 --- a/OGS/rf.cpp +++ b/OGS/rf.cpp @@ -227,16 +227,6 @@ int main(int argc, char* argv[]) std::cout << BuildInfo::OGS_VERSION << "\n"; exit(0); } - if (anArg == "--model-root" || anArg == "-m") - { - if (i + 1 >= arg_strings.size()) - { - std::cerr << "Error: Parameter " << anArg << " needs an additional argument" << std::endl; - std::exit(EXIT_FAILURE); - } -// modelRoot = std::string(arg_strings[++i]); - continue; - } if (anArg == "--output-directory") { if (i + 1 >= arg_strings.size())