From d831aa99a3c282349260010ea22502aad7acf697 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 5 Jun 2024 15:58:15 +1000 Subject: [PATCH 1/8] Adding auto timestep to lammps integration --- src/lammps.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/lammps.cpp b/src/lammps.cpp index 987ec796..31f10996 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -121,6 +121,7 @@ int main(int argc, char* argv[]) { exit(1); } if (RFI.hasVar("output")) { + fprintf(fp, "variable timestep equal %.15lg\n", RFI.auto_timestep); fprintf(fp, "variable output string %s\n", RFI.getVar("output").c_str()); } fprintf(fp, "%s\n", RFI.getVar("content").c_str()); From 51f3ea19b17ddf649f6987ad502500af0d16988e Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 5 Jun 2024 15:58:48 +1000 Subject: [PATCH 2/8] Adding TCLB/LIGGGHTS example --- example/data/pipe2x1.stl | Bin 0 -> 3084 bytes example/particle/3d/cylinder.xml | 85 +++++++++++++++++++++++++++++++ 2 files changed, 85 insertions(+) create mode 100644 example/data/pipe2x1.stl create mode 100644 example/particle/3d/cylinder.xml diff --git a/example/data/pipe2x1.stl b/example/data/pipe2x1.stl new file mode 100644 index 0000000000000000000000000000000000000000..ab8f2b40708a59b3259bb37f9d41cbb830239705 GIT binary patch literal 3084 zcma)8F>4e-6rRREprArjuGz{;n>1#2Y^)9x4CyQmJdVX9hvtOj3R!=Et)#NDvrN6a zwV=d8)W#Gk0@h+9#`xa5_jbOQYg{<^_P+1iZ)Rua?c1t$o{y&0Vl=O2+k?&7Z2F{{ z4<_U7@w3r<{CMzicz1BGnvTPNUUyV?y8oj5b@59>w*~!5{W}3Jjwiq4$CJ;jp$i-d z7*r*a*^>y)c6Z8`>u*|uPKEC4M?d2G@&Rk;0!IR-0F2q>_;&E9{Bw3D8K?qWG;4?P z>fUS3KnF0Gm*9WOXksjj;mwX9J_}VvF?<)B*lGqi5-_MjH$C{u-sQ4cyDu5wKsTxC z`1RnsW`HAkrdOSFAC{luwWEqN=*lX- zg-cr-r~7|30~`q%z%Dz~JbFr%3{(N;j3#G*gEM*tu*(iLjjFYIoKy8~np2zK6*zdj zJp1?sP9P#ynep zb?8nkc^>EhYxV2j>#&1S9w&{ zGe`zF&}|)uW@yhK^Pt_w9Sbw4LN`&`$*{wz&`rCY`zy@=2cI=#taD-&?_JBLTkXhX zB;asvI`iBoE;F117?Z>uXBhT)CY~biX literal 0 HcmV?d00001 diff --git a/example/particle/3d/cylinder.xml b/example/particle/3d/cylinder.xml new file mode 100644 index 00000000..dbbebe77 --- /dev/null +++ b/example/particle/3d/cylinder.xml @@ -0,0 +1,85 @@ + + + + + + + + + + + + + + + + + + + + + +units cgs +boundary p f f +newton off # required off for tangential history +atom_style sphere +atom_modify map array +atom_modify sort 1 0.4 +communicate single vel yes +processors * 1 1 + +neighbor 0.006 bin # ensure skin distance + rp_lrg + rp_sml > dp_lrg +neigh_modify delay 0 + +# Declare domain +region domain block 0 0.18 -0.045 0.045 -0.045 0.045 +create_box 1 domain + +# Specify particle groups +group particle_group type 1 + +# Define region for particle insertion +region pack cylinder x 0 0 0.0375 0 0.12 + +# Insert particles +fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant 0.0075000000 +fix dist particle_group particledistribution/discrete 18143 1 part_1 1 +fix ins particle_group insert/pack seed 100003 distributiontemplate dist maxattempt 500 insert_every once overlapcheck yes all_in yes region pack volumefraction_region 0.050000 check_dist_from_subdomain_border no +run 1 + +# Specify particle groups +group particle_group type 1 + +# Define material properties (from which kn kt etc. are calculated for hertz interactions) +soft_particles yes +fix m1 all property/global youngsModulus peratomtype 3000.000000 # defines kn, kt, gamma_n, gamma_t +fix m2 all property/global poissonsRatio peratomtype 0.5 # defines kn, kt, gamma_n, gamma_t +fix m3 all property/global coefficientRestitution peratomtypepair 1 0.8 # defines damping, must be >0.05 +fix m4 all property/global coefficientFriction peratomtypepair 1 0.5 # defines friction + +fix frac_wall all mesh/surface file example/data/pipe2x1.stl type 1 scale 0.09 move 0 0 0 + +# Define physics for particle interactions +pair_style gran model hertz tangential history # 'tangential off' sets Ft=0; 'tangential no_history' incorporates damping to Ft, sets kt=0; 'tangential history' incorporate kt and damping into Ft +pair_coeff * * + +fix granwalls all wall/gran model hertz tangential history mesh n_meshes 1 meshes frac_wall + +# Apply integration +fix integr particle_group nve/sphere + +# Couple to TCLB +fix tclb all external pf/callback 1 1 + +dump vtk_dump all atom/vtk 1000 ${output}_part_*.vtu + +timestep ${timestep} + +run 50000 + + + + + + + From 1335c00ec4c0270a37badf4b2badd566eb92de21 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 5 Jun 2024 16:03:43 +1000 Subject: [PATCH 3/8] Correcting error in writing file in multiple threads --- src/lammps.cpp | 24 +++++++++++++----------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/src/lammps.cpp b/src/lammps.cpp index 31f10996..13daf1f3 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -113,19 +113,21 @@ int main(int argc, char* argv[]) { MPI_Abort(MPI_COMM_WORLD, 1); exit(1); } - FILE* fp = NULL; - fp = fopen(infile.c_str(), "w"); - if (fp == NULL) { - printf("ERROR: failed to write to %s\n",infile.c_str()); - MPI_Abort(MPI_COMM_WORLD, 1); - exit(1); - } - if (RFI.hasVar("output")) { + if (MPMD.local_rank == 0) { + FILE* fp = NULL; + fp = fopen(infile.c_str(), "w"); + if (fp == NULL) { + printf("ERROR: failed to write to %s\n",infile.c_str()); + MPI_Abort(MPI_COMM_WORLD, 1); + exit(1); + } fprintf(fp, "variable timestep equal %.15lg\n", RFI.auto_timestep); - fprintf(fp, "variable output string %s\n", RFI.getVar("output").c_str()); + if (RFI.hasVar("output")) { + fprintf(fp, "variable output string %s\n", RFI.getVar("output").c_str()); + } + fprintf(fp, "%s\n", RFI.getVar("content").c_str()); + fclose(fp); } - fprintf(fp, "%s\n", RFI.getVar("content").c_str()); - fclose(fp); } else { printf("ERROR: No configuration provided (either xml file or content from force calculator\n"); MPI_Abort(MPI_COMM_WORLD,1); From 4992601620c6f95ca0d8764b4fb7d4aa4c937f20 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 5 Jun 2024 16:52:55 +1000 Subject: [PATCH 4/8] First implementation of velocity ease-in for simplepart --- src/simplepart.cpp | 41 ++++++++++++++++++++++++++++++++++++----- 1 file changed, 36 insertions(+), 5 deletions(-) diff --git a/src/simplepart.cpp b/src/simplepart.cpp index e5f5b582..cedcecba 100644 --- a/src/simplepart.cpp +++ b/src/simplepart.cpp @@ -6,16 +6,20 @@ #include const double twopi = 8*atan(1.0); +const double pi = 4*atan(1.0); struct Particle { double x[3]; double r; double m; double v[3]; + double v0[3]; double f[3]; double favg[3]; double omega[3]; + double omega0[3]; double torque[3]; + double ease_in_time; size_t n; bool logging; Particle() { @@ -31,6 +35,7 @@ struct Particle { m = 0; r = 0; logging = false; + ease_in_time = 0; } }; @@ -180,15 +185,17 @@ int main(int argc, char *argv[]) { if (attr_name.vector == "") { p.x[attr_name.d] = attr.as_double(); } else if (attr_name.vector == "v") { - p.v[attr_name.d] = attr.as_double(); + p.v0[attr_name.d] = attr.as_double(); } else if (attr_name.vector == "omega") { - p.omega[attr_name.d] = attr.as_double(); + p.omega0[attr_name.d] = attr.as_double(); } else if (attr_name == "r") { p.r = attr.as_double(); } else if (attr_name == "m") { p.m = attr.as_double(); } else if (attr_name == "log") { p.logging = attr.as_bool(); + } else if (attr_name == "ease-in") { + p.ease_in_time = attr.as_double(); } else { ERROR("Unknown atribute '%s' in '%s'", attr.name(), node.name()); return -1; @@ -198,6 +205,17 @@ int main(int argc, char *argv[]) { ERROR("Specify the radius with 'r' attribute"); return -1; } + if (p.ease_in_time > 0) { + for (int i=0;i<3;i++) { + p.v[i] = 0; + p.omega[i] = 0; + } + } else { + for (int i=0;i<3;i++) { + p.v[i] = p.v0[i]; + p.omega[i] = p.omega0[i]; + } + } p.n = particles.size(); particles.push_back(p); } else if (node_name == "Periodic") { @@ -387,11 +405,24 @@ int main(int argc, char *argv[]) { fprintf(logging_f, "\n"); } for (Particles::iterator p = particles.begin(); p != particles.end(); p++) { + double t = dt * iter; for (int i=0; i<3; i++) p->favg[i] = p->favg[i] + p->f[i]; - if (p->m > 0.0) { - for (int i=0; i<3; i++) p->v[i] = p->v[i] + p->f[i] / p->m * dt; + if (p->ease_in_time > 0) { + if (p->ease_in_time > t) { + double fac = (1-cos(pi * t / p->ease_in_time))*0.5; + for (int i=0; i<3; i++) p->v[i] = p->v0[i] * fac; + for (int i=0; i<3; i++) p->omega[i] = p->omega0[i] * fac; + } else { + for (int i=0; i<3; i++) p->v[i] = p->v0[i]; + for (int i=0; i<3; i++) p->omega[i] = p->omega0[i]; + p->ease_in_time = 0; + } + } else { + if (p->m > 0.0) { + for (int i=0; i<3; i++) p->v[i] = p->v[i] + p->f[i] / p->m * dt; + } + for (int i=0; i<3; i++) p->v[i] = p->v[i] + acc_vec[i] * cos(twopi * t * acc_freq); } - for (int i=0; i<3; i++) p->v[i] = p->v[i] + acc_vec[i] * cos(twopi * dt * iter * acc_freq); for (int i=0; i<3; i++) p->x[i] = p->x[i] + p->v[i] * dt; } iter++; From 62c29d47f6157097557135596db5439a7f27ea18 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 19 Jun 2024 14:09:40 +1000 Subject: [PATCH 5/8] Adding first version of shear cell example --- example/data/plane1x1.stl | Bin 0 -> 184 bytes example/particle/3d/shear1.xml | 82 +++++++++++++++++++++++++++++++++ 2 files changed, 82 insertions(+) create mode 100644 example/data/plane1x1.stl create mode 100644 example/particle/3d/shear1.xml diff --git a/example/data/plane1x1.stl b/example/data/plane1x1.stl new file mode 100644 index 0000000000000000000000000000000000000000..d528b01d46f0e7c3255e9f31b88a0d42692e71af GIT binary patch literal 184 zcmWH`EG|vV$*f8&$;{7F2+7aS$<8cMNKeg6ElMm&O;HH;aa3^2%t>V+5-_2tYp_S* WGT`AiKvaTMz%WcNOc#s>i30$~aTQ_! literal 0 HcmV?d00001 diff --git a/example/particle/3d/shear1.xml b/example/particle/3d/shear1.xml new file mode 100644 index 00000000..15ea7b3a --- /dev/null +++ b/example/particle/3d/shear1.xml @@ -0,0 +1,82 @@ + + + + + + + + + + + + + + + + + units cgs + boundary p f f + newton off # required off for tangential history + atom_style sphere + atom_modify map array + atom_modify sort 1 0.4 + communicate single vel yes + processors * 1 1 + + neighbor 0.006 bin # ensure skin distance + rp_lrg + rp_sml > dp_lrg + neigh_modify delay 0 + + # Declare domain + region domain block -1 1 0 1 -1 1 + create_box 1 domain + + # Specify particle groups + group particle_group type 1 + + # Define region for particle insertion + region pack block -1 1 0 1 -1 1 + + # Insert particles + fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant 0.1 + fix dist particle_group particledistribution/discrete 18143 1 part_1 1 + fix ins particle_group insert/pack seed 100003 distributiontemplate dist maxattempt 500 insert_every once overlapcheck yes all_in yes region pack volumefraction_region 0.050000 check_dist_from_subdomain_border no + run 1 + + # Specify particle groups + group particle_group type 1 + + # Define material properties (from which kn kt etc. are calculated for hertz interactions) + soft_particles yes + fix m1 all property/global youngsModulus peratomtype 3000.000000 # defines kn, kt, gamma_n, gamma_t + fix m2 all property/global poissonsRatio peratomtype 0.5 # defines kn, kt, gamma_n, gamma_t + fix m3 all property/global coefficientRestitution peratomtypepair 1 0.8 # defines damping, must be >0.05 + fix m4 all property/global coefficientFriction peratomtypepair 1 0.5 # defines friction + + fix top_wall all mesh/surface file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 1 -1 surface_vel 1 0 0 + fix bottom_wall all mesh/surface file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 0 -1 + + # Define physics for particle interactions + pair_style gran model hertz tangential history # 'tangential off' sets Ft=0; 'tangential no_history' incorporates damping to Ft, sets kt=0; 'tangential history' incorporate kt and damping into Ft + pair_coeff * * + + fix granwalls all wall/gran model hertz tangential history mesh n_meshes 2 meshes top_wall bottom_wall + + # Apply integration + fix integr particle_group nve/sphere + + # Couple to TCLB + fix tclb all external pf/callback 1 1 + + dump vtk_dump all atom/vtk 1000 ${output}_part_*.vtu + + timestep ${timestep} + + run 50000 + + + + + + + + From 26c5da52d09291824f45abc49e2d90478edafb35 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Wed, 19 Jun 2024 14:21:04 +1000 Subject: [PATCH 6/8] Adding particle force measurement --- example/particle/3d/shear1.xml | 23 +++++++++++++++++------ 1 file changed, 17 insertions(+), 6 deletions(-) diff --git a/example/particle/3d/shear1.xml b/example/particle/3d/shear1.xml index 15ea7b3a..1278dd94 100644 --- a/example/particle/3d/shear1.xml +++ b/example/particle/3d/shear1.xml @@ -8,11 +8,11 @@ - + - + units cgs boundary p f f @@ -39,7 +39,7 @@ # Insert particles fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant 0.1 fix dist particle_group particledistribution/discrete 18143 1 part_1 1 - fix ins particle_group insert/pack seed 100003 distributiontemplate dist maxattempt 500 insert_every once overlapcheck yes all_in yes region pack volumefraction_region 0.050000 check_dist_from_subdomain_border no + fix ins particle_group insert/pack seed 100003 distributiontemplate dist maxattempt 500 insert_every once overlapcheck yes all_in yes region pack volumefraction_region 0.30000 check_dist_from_subdomain_border no run 1 # Specify particle groups @@ -52,14 +52,14 @@ fix m3 all property/global coefficientRestitution peratomtypepair 1 0.8 # defines damping, must be >0.05 fix m4 all property/global coefficientFriction peratomtypepair 1 0.5 # defines friction - fix top_wall all mesh/surface file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 1 -1 surface_vel 1 0 0 - fix bottom_wall all mesh/surface file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 0 -1 + fix topwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 1 -1 surface_vel 1 0 0 + fix bottomwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 0 -1 # Define physics for particle interactions pair_style gran model hertz tangential history # 'tangential off' sets Ft=0; 'tangential no_history' incorporates damping to Ft, sets kt=0; 'tangential history' incorporate kt and damping into Ft pair_coeff * * - fix granwalls all wall/gran model hertz tangential history mesh n_meshes 2 meshes top_wall bottom_wall + fix granwalls all wall/gran model hertz tangential history mesh n_meshes 2 meshes topwall bottomwall # Apply integration fix integr particle_group nve/sphere @@ -67,8 +67,19 @@ # Couple to TCLB fix tclb all external pf/callback 1 1 + variable time equal step*dt + variable tfx equal f_topwall[1] + variable tfy equal f_topwall[2] + variable tfz equal f_topwall[3] + variable bfx equal f_bottomwall[1] + variable bfy equal f_bottomwall[2] + variable bfz equal f_bottomwall[3] + dump forces all mesh/vtk 1000 ${output}_wall_*.vtk output interpolate id stress stresscomponents + fix forceslog all print 100 "${time},${tfx},${tfy},${tfz},${bfx},${bfy},${bfz}" file ${output}_forces.csv title "t,tFx,tFy,tFz,bFx,bFy,bFz" screen no + dump vtk_dump all atom/vtk 1000 ${output}_part_*.vtu + timestep ${timestep} run 50000 From a3ccf75485512387c32d1131e03b6227eecc10f9 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Mon, 24 Jun 2024 12:11:38 +1000 Subject: [PATCH 7/8] Passing general attributes through RFI --- example/particle/3d/shear1.xml | 4 +- src/Handlers/acRemoteForceInterface.cpp | 55 ++++++++++++++----------- src/RemoteForceInterface.h | 5 +++ src/lammps.cpp | 21 ++++++++-- 4 files changed, 55 insertions(+), 30 deletions(-) diff --git a/example/particle/3d/shear1.xml b/example/particle/3d/shear1.xml index 1278dd94..000e9147 100644 --- a/example/particle/3d/shear1.xml +++ b/example/particle/3d/shear1.xml @@ -13,7 +13,7 @@ - + units cgs boundary p f f newton off # required off for tangential history @@ -37,7 +37,7 @@ region pack block -1 1 0 1 -1 1 # Insert particles - fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant 0.1 + fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant ${radius} fix dist particle_group particledistribution/discrete 18143 1 part_1 1 fix ins particle_group insert/pack seed 100003 distributiontemplate dist maxattempt 500 insert_every once overlapcheck yes all_in yes region pack volumefraction_region 0.30000 check_dist_from_subdomain_border no run 1 diff --git a/src/Handlers/acRemoteForceInterface.cpp b/src/Handlers/acRemoteForceInterface.cpp index c1198245..ac644a59 100644 --- a/src/Handlers/acRemoteForceInterface.cpp +++ b/src/Handlers/acRemoteForceInterface.cpp @@ -15,7 +15,6 @@ int acRemoteForceInterface::Init () { int acRemoteForceInterface::ConnectRemoteForceInterface(std::string integrator_) { output("Connecting RFI to %s\n",integrator_.c_str()); - pugi::xml_attribute attr; double units[3]; units[0] = solver->units.alt("1m"); units[1] = solver->units.alt("1s"); @@ -25,7 +24,6 @@ int acRemoteForceInterface::ConnectRemoteForceInterface(std::string integrator_) solver->lattice->RFI.CanCopeWithUnits(false); solver->lattice->RFI.setVar("output", solver->info.outpath); - std::string element_content; int node_children = 0; @@ -44,24 +42,40 @@ int acRemoteForceInterface::ConnectRemoteForceInterface(std::string integrator_) } } if (node_children > 0) solver->lattice->RFI.setVar("content", element_content); + bool stats = false; std::string stats_prefix = solver->info.outpath; stats_prefix = stats_prefix + "_RFI"; int stats_iter = 200; - - attr = node.attribute("stats"); - if (attr) stats = attr.as_bool(); - attr = node.attribute("stats_iter"); - if (attr) { - stats_iter = solver->units.alt(attr.value()); - stats = true; - } - attr = node.attribute("stats_prefix"); - if (attr) { - stats_prefix = attr.value(); - stats = true; + bool use_box = true; + + + for (pugi::xml_attribute attr = node.first_attribute(); attr; attr = attr.next_attribute()) { + std::string attr_name = attr.name(); + if (attr_name == "integrator") { + // ignore + } else if (attr_name == "stats") { + stats = attr.as_bool(); + } else if (attr_name == "stats_iter") { + stats_iter = solver->units.alt(attr.value()); + stats = true; + } else if (attr_name == "stats_prefix") { + stats_prefix = attr.value(); + stats = true; + } else if (attr_name == "use_box") { + use_box = attr.as_bool(); + } else if (attr_name == "omega") { + solver->lattice->RFI_omega = attr.as_bool(); + } else if (attr_name == "torque") { + solver->lattice->RFI_torque = attr.as_bool(); + } else { + double val = solver->units.alt(attr.value()); + char str[STRING_LEN]; + sprintf(str, "%.15lg", val); + solver->lattice->RFI.setVar(attr.name(), str); + } } - + if (stats) { output("Asking for stats on RFI ( %s every %d it)\n", stats_prefix.c_str(), stats_iter); solver->lattice->RFI.enableStats(stats_prefix.c_str(), stats_iter); @@ -74,10 +88,6 @@ int acRemoteForceInterface::ConnectRemoteForceInterface(std::string integrator_) } integrator = integrator_; - bool use_box = true; - attr = node.attribute("use_box"); - if (attr) use_box = attr.as_bool(); - if (use_box) { lbRegion reg = solver->lattice->region; double px = solver->lattice->px; @@ -92,15 +102,10 @@ int acRemoteForceInterface::ConnectRemoteForceInterface(std::string integrator_) pz + reg.dz + reg.nz + PART_MAR_BOX); } - attr = node.attribute("omega"); - if (attr) solver->lattice->RFI_omega = attr.as_bool(); - attr = node.attribute("torque"); - if (attr) solver->lattice->RFI_torque = attr.as_bool(); - MPI_Barrier(MPMD.local); solver->lattice->RFI.Connect(MPMD.work,inter.work); - return 0; + return 0; } diff --git a/src/RemoteForceInterface.h b/src/RemoteForceInterface.h index 1dc10912..9d561af5 100644 --- a/src/RemoteForceInterface.h +++ b/src/RemoteForceInterface.h @@ -151,6 +151,11 @@ class RemoteForceInterface { void setUnits(rfi_real_t meter, rfi_real_t second, rfi_real_t kilogram); void setVar(const vars_name_t& name, const vars_value_t& value); bool hasVar(const vars_name_t& name) { return vars.find(name) != vars.end(); }; + std::vector listVars() { + std::vector names; + for (auto it : vars) names.push_back(it.first); + return names; + } const vars_value_t& getVar(const vars_name_t& name) { return vars[name]; }; inline rfi_real_t& RawData(size_t i, int j) { if (STORAGE == ArrayOfStructures) { diff --git a/src/lammps.cpp b/src/lammps.cpp index 13daf1f3..d57379ff 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -121,11 +121,26 @@ int main(int argc, char* argv[]) { MPI_Abort(MPI_COMM_WORLD, 1); exit(1); } - fprintf(fp, "variable timestep equal %.15lg\n", RFI.auto_timestep); - if (RFI.hasVar("output")) { - fprintf(fp, "variable output string %s\n", RFI.getVar("output").c_str()); + const std::vector var_names = RFI.listVars(); + for (const auto& v : var_names) { + if (v == "content") continue; + auto& value = RFI.getVar(v); + bool is_numeric = false; + if (v != "output") { + double val; + int ret, len; + ret = sscanf(value.c_str(),"%lf%n", &val, &len); + if ((ret > 0) && (len == value.size())) { + is_numeric = true; + fprintf(fp, "variable %s equal %.15lg\n", v.c_str(), val); + } + } + if (!is_numeric) fprintf(fp, "variable %s string %s\n", v.c_str(), value.c_str()); } + fprintf(fp, "variable timestep equal %.15lg\n", RFI.auto_timestep); + fprintf(fp, "\n"); fprintf(fp, "%s\n", RFI.getVar("content").c_str()); + fprintf(fp, "\n"); fclose(fp); } } else { From 3c1f54cc61b5b907642348952bd202f32a7ec0f9 Mon Sep 17 00:00:00 2001 From: L Laniewski-Wollk Date: Mon, 24 Jun 2024 15:12:55 +1000 Subject: [PATCH 8/8] Updating the shear1 example --- example/particle/3d/shear1.xml | 29 +++++++++++++++-------------- src/lammps.cpp | 4 ++-- 2 files changed, 17 insertions(+), 16 deletions(-) diff --git a/example/particle/3d/shear1.xml b/example/particle/3d/shear1.xml index 000e9147..1bc1e339 100644 --- a/example/particle/3d/shear1.xml +++ b/example/particle/3d/shear1.xml @@ -2,8 +2,8 @@ - - + + @@ -13,7 +13,7 @@ - + units cgs boundary p f f newton off # required off for tangential history @@ -27,14 +27,15 @@ neigh_modify delay 0 # Declare domain - region domain block -1 1 0 1 -1 1 + variable len2 equal ${length}*0.5 + region domain block -${len2} ${len2} 0 ${height} -${len2} ${len2} create_box 1 domain # Specify particle groups group particle_group type 1 # Define region for particle insertion - region pack block -1 1 0 1 -1 1 + region pack block -${len2} ${len2} 0 ${height} -${len2} ${len2} # Insert particles fix part_1 particle_group particletemplate/sphere 17891 atom_type 1 density constant 1.0 radius constant ${radius} @@ -47,13 +48,13 @@ # Define material properties (from which kn kt etc. are calculated for hertz interactions) soft_particles yes - fix m1 all property/global youngsModulus peratomtype 3000.000000 # defines kn, kt, gamma_n, gamma_t + fix m1 all property/global youngsModulus peratomtype 30000.000000 # defines kn, kt, gamma_n, gamma_t fix m2 all property/global poissonsRatio peratomtype 0.5 # defines kn, kt, gamma_n, gamma_t fix m3 all property/global coefficientRestitution peratomtypepair 1 0.8 # defines damping, must be >0.05 fix m4 all property/global coefficientFriction peratomtypepair 1 0.5 # defines friction - fix topwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 1 -1 surface_vel 1 0 0 - fix bottomwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale 2 rotate axis 1 0 0 angle 90 move -1 0 -1 + fix topwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale ${length} rotate axis 1 0 0 angle 90 move -${len2} ${height} -${len2} surface_vel 1 0 0 + fix bottomwall all mesh/surface/stress file example/data/plane1x1.stl type 1 scale ${length} rotate axis 1 0 0 angle 90 move -${len2} 0 -${len2} # Define physics for particle interactions pair_style gran model hertz tangential history # 'tangential off' sets Ft=0; 'tangential no_history' incorporates damping to Ft, sets kt=0; 'tangential history' incorporate kt and damping into Ft @@ -74,20 +75,20 @@ variable bfx equal f_bottomwall[1] variable bfy equal f_bottomwall[2] variable bfz equal f_bottomwall[3] - dump forces all mesh/vtk 1000 ${output}_wall_*.vtk output interpolate id stress stresscomponents - fix forceslog all print 100 "${time},${tfx},${tfy},${tfz},${bfx},${bfy},${bfz}" file ${output}_forces.csv title "t,tFx,tFy,tFz,bFx,bFy,bFz" screen no + dump forces all mesh/vtk ${vtk_it} ${output}_wall_*.vtk output interpolate id stress stresscomponents + fix forceslog all print ${log_it} "${time},${tfx},${tfy},${tfz},${bfx},${bfy},${bfz}" file ${output}_forces.csv title "t,tFx,tFy,tFz,bFx,bFy,bFz" screen no - dump vtk_dump all atom/vtk 1000 ${output}_part_*.vtu + dump vtk_dump all atom/vtk ${vtk_it} ${output}_part_*.vtu timestep ${timestep} - run 50000 + run ${iterations} - + - + diff --git a/src/lammps.cpp b/src/lammps.cpp index d57379ff..44d04939 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -132,12 +132,12 @@ int main(int argc, char* argv[]) { ret = sscanf(value.c_str(),"%lf%n", &val, &len); if ((ret > 0) && (len == value.size())) { is_numeric = true; - fprintf(fp, "variable %s equal %.15lg\n", v.c_str(), val); + fprintf(fp, "variable %s equal %.13lg\n", v.c_str(), val); } } if (!is_numeric) fprintf(fp, "variable %s string %s\n", v.c_str(), value.c_str()); } - fprintf(fp, "variable timestep equal %.15lg\n", RFI.auto_timestep); + fprintf(fp, "variable timestep equal %.13lg\n", RFI.auto_timestep); fprintf(fp, "\n"); fprintf(fp, "%s\n", RFI.getVar("content").c_str()); fprintf(fp, "\n");