Added restart functionality for bonds

Testfile: examples/bondspackage/beam/read_restart.bonds

examples/bondspackage/beam/beam.bonds

@@ -144,3 +144,5 @@ fix 	tipF tipball addforce  0 0 -$F    #1kip=4.4482216kN
 variable t equal 1/dt 	#run 1second
 run $t
 
+write_restart beam.restart
+

examples/bondspackage/beam/read_restart.bonds

@@ -0,0 +1,138 @@
+#bond
+
+######## simulation setup #########
+variable dt equal 0.00001
+variable runstep equal 10/${dt}
+variable filestep equal 0.01/${dt}
+##################################
+
+variable rp equal 10*0.0254   #particle radius
+variable d  equal 2*${rp}
+variable sc equal $d-0.000000001
+variable c  equal ${d}+0.1
+variable rb equal 0.8*${rp}
+variable L  equal 200*0.0254+${rp}	    #INCH !
+
+atom_style	hybrid granular bond/gran n_bondtypes 1 bonds_per_atom 6
+#atom_style	granular
+atom_modify	map array
+boundary	f f f
+newton		off
+
+#processors  2 1 1
+
+communicate	single vel yes
+
+units		si
+
+region		reg block -$L $L -1 1 -$L 1 units box
+#create_box	2 reg #bonds 1 6 			# {num_atom_types} {region} bonds moved into atom_style for LIGGGHTS>=3.X.X
+
+#bond
+bond_style 	gran
+
+#read_restart can only work correctly AFTER bond_style is set
+read_restart beam.restart
+
+neighbor	1 bin
+neigh_modify	delay 0
+
+#New pair style
+pair_style 	gran model hertz tangential history 
+pair_coeff	* *
+
+#for stressbreak
+variable Y  equal 6894.75729*1000*1e4
+
+variable p equal 0.3
+variable G equal 3759.4*6894.75729*1000 #$Y/(2*(1+$p))
+print $Y
+print $G
+
+#double kn = E/lenght; 				 //ZUG - Druck
+#double ks = 12*I*E/(A*lenght*lenght*lenght);    //A*G/lenght; //Schub
+
+variable I equal 0.25*PI*${rb}*${rb}*${rb}*${rb}
+variable A equal PI*${rb}*${rb}
+
+variable kn equal 2*$Y/$d 		 #grain behaivior ?
+variable ks equal 12*$I*$Y/($A*$d*$d*$d) #${kn}/2.5
+
+print $I
+print $A
+print ${kn}
+print ${ks}
+bond_coeff 	1 ${rb} ${kn} ${ks} 1 1e10 1e10
+#bond_coeff 	1 ${rb} $Y $G 1 1e10 1e10 #N rb Sn_bond St_bond TYPE_OF_BOND r_break|sigma_break tau_break (temp)
+
+
+#Material properties required for new pair styles
+
+fix 		m1 all property/global youngsModulus peratomtype $Y $Y
+fix 		m2 all property/global poissonsRatio peratomtype $p $p
+fix 		m3 all property/global coefficientRestitution peratomtypepair 2 0.2 0.2 0.2 0.2
+fix 		m4 all property/global coefficientFriction peratomtypepair 2 0.3 0.3 0.3 0.3
+
+mass 		1 1.0 #dummy
+mass 		2 1.0 #dummy
+
+variable lastone equal atoms
+
+group		freezed id == 1
+group		tipball id == ${lastone}
+
+run 1
+fix f1	freezed freeze	#festlager
+
+fix 		bondcr all bond/create/gran 1 1 1 0.51 1 6 #every itype jtype cutoff btype newperts
+
+timestep	${dt}
+
+#fix		gravi all gravity 9.81 vector 0.0 0.0 -1.0
+
+#fix		zwalls all wall/gran/hertz/history primitive type 1 zplane 0.0 
+
+#region and insertion
+
+group		nve_group region reg
+
+#apply nve integration to all particles that are inserted as single particles
+fix		integr nve_group nve/sphere
+
+#output settings, include total thermal energyg
+thermo_style	custom step atoms #numbond
+thermo		100
+thermo_modify	lost ignore norm no
+compute_modify	thermo_temp dynamic yes
+
+#insert the first particles so that dump is not empty
+shell rm -rf post
+shell mkdir post
+dump		dmp 	all custom  ${filestep} post/dump*.liggghts id type x y z vx vy vz fx fy fz omegax omegay omegaz radius
+
+compute b1 all property/local batom1x batom1y batom1z batom2x batom2y batom2z batom1 batom2 btype
+dump 	bnd all local ${filestep} post/bonds*.bond c_b1[1] c_b1[2] c_b1[3] c_b1[4] c_b1[5] c_b1[6] c_b1[7] c_b1[8] c_b1[9]  #x1 y1 z1 x2 y2 z2 id1 id2 typeofbond BSD
+
+#hole koordinaten von letztem partikel
+variable time equal step*dt
+variable px equal x[${lastone}]
+variable py equal y[${lastone}]
+variable pz equal z[${lastone}]
+
+fix outfile all print ${filestep} "${time};${px};${py};${pz}" file beam.csv screen no title "t;x;y;z"
+
+#insert particles
+run		1		#run one step
+fix_modify	bondcr every 0  #do not create new bonds after this line
+variable KIP equal 4448.2216 #N
+variable inch equal 0.0254   #mm
+
+variable M equal 100*${KIP}/${inch}
+variable F equal 1*${KIP}		#N
+fix 	tipF tipball addforce  0 0 +$F    #1kip=4.4482216kN
+#fix 	tipM tipball addtorque 0 $M  0	  #nicht moeglich 
+
+#quit
+variable t equal 1/dt 	#run 1second
+run $t
+

src/atom_vec_bond_gran.cpp

@@ -66,7 +66,7 @@ AtomVecBondGran::AtomVecBondGran(LAMMPS *lmp) : AtomVec(lmp)
 
 void AtomVecBondGran::settings(int narg, char **arg)
 {
-
+  if (narg == 0) return;		//in case of restart no arguments are given, instead nbondtypes and bond_per_atom are defined by read_restart_settings
   if (narg != 4) error->all(FLERR,"Invalid atom_style bond/gran command,expecting exactly 4 arguments");
 
   if(strcmp(arg[0],"n_bondtypes"))
@@ -107,7 +107,7 @@ void AtomVecBondGran::grow(int n)
   else nmax = n;
   atom->nmax = nmax;
 
-    tag = memory->grow(atom->tag,nmax,"atom:tag");
+  tag = memory->grow(atom->tag,nmax,"atom:tag");
   type = memory->grow(atom->type,nmax,"atom:type");
   mask = memory->grow(atom->mask,nmax,"atom:mask");
   image = memory->grow(atom->image,nmax,"atom:image");
@@ -121,6 +121,12 @@ void AtomVecBondGran::grow(int n)
   bond_type = memory->grow(atom->bond_type,nmax,atom->bond_per_atom,"atom:bond_type");
   bond_atom = memory->grow(atom->bond_atom,nmax,atom->bond_per_atom,"atom:bond_atom");
 
+  if(0 == atom->bond_per_atom)
+    error->all(FLERR,"Bonded particles need bond_per_atom > 0");
+
+  if(atom->n_bondhist < 0)
+	  error->all(FLERR,"atom->n_bondhist < 0 suggests that 'bond_style gran' has not been called before 'read_restart' command! Please check that.");
+
   if(atom->n_bondhist)
   {
      bond_hist = atom->bond_hist =
@@ -657,7 +663,7 @@ int AtomVecBondGran::size_restart()
   {
     n += 13 + 2*num_bond[i];
 
-    if(atom->n_bondhist) n += 1/*num_bondhist*/ + num_bond[i] * num_bondhist/*bond_hist*/;
+    if(atom->n_bondhist) n += 1/*num_bondhist*/ + num_bond[i] * atom->n_bondhist/*bond_hist*/; //CR 26.01.2015
   }
 
   if (atom->nextra_restart)
@@ -753,9 +759,10 @@ int AtomVecBondGran::unpack_restart(double *buf)
   {
       if(atom->n_bondhist != static_cast<int>(buf[m++]))
           error->all(FLERR,"Íncompatibel restart file: file was created using a bond model with a different number of history values");
-      for (k = 0; k < num_bond[nlocal]; k++)
+     printf("num_bond[nlocal]=%d\n",num_bond[nlocal]); 
+     for (k = 0; k < num_bond[nlocal]; k++)
          for (l = 0; l < atom->n_bondhist; l++)
-            bond_hist[nlocal][k][l] = buf[m++];
+            atom->bond_hist[nlocal][k][l] = buf[m++];
   }
 
   double **extra = atom->extra;
@@ -868,6 +875,22 @@ void AtomVecBondGran::write_data(FILE *fp, int n, double **buf)
   error->all(FLERR,"Add usefull code here");
 }
 
+void AtomVecBondGran::write_restart_settings(FILE *fp)
+{
+  fwrite(&atom->nbondtypes,sizeof(int),1,fp);
+  fwrite(&atom->bond_per_atom,sizeof(int),1,fp);
+}
+
+void AtomVecBondGran::read_restart_settings(FILE *fp)
+{
+  if (comm->me == 0) {
+    fread(&atom->nbondtypes,sizeof(int),1,fp);
+    fread(&atom->bond_per_atom,sizeof(int),1,fp);
+  }
+  MPI_Bcast(&atom->nbondtypes,1,MPI_INT,0,world);
+  MPI_Bcast(&atom->bond_per_atom,1,MPI_INT,0,world);
+}
+
 /* ----------------------------------------------------------------------
    return # of bytes of allocated memory
 ------------------------------------------------------------------------- */

src/atom_vec_bond_gran.h

@@ -65,6 +65,8 @@ class AtomVecBondGran : public AtomVec {
   void pack_data(double **);
   void pack_data(double **buf,int tag_offset); 
   void write_data(FILE *, int, double **);
+  void write_restart_settings(FILE *);
+  void read_restart_settings(FILE *);
 
  private:
   int *tag,*type,*mask,*image;

src/fix_bond_propagate_gran.cpp

@@ -211,7 +211,7 @@ inline void FixBondPropagateGran::remove_bond(int ilocal,int ibond, int bondnumb
 
 void FixBondPropagateGran::write_restart(FILE *fp)
 {
-  error->warning(FLERR,"Restart functionality not yet tested for granular bonds...");
+  //error->warning(FLERR,"Restart functionality not yet tested for granular bonds...");
 
   //NP write a dummy value
   int n = 0;
@@ -240,5 +240,5 @@ void FixBondPropagateGran::restart(char *buf)
 
   double dummy = static_cast<int> (list[n++]);
 
-  error->warning(FLERR,"Restart functionality not yet tested for granular bonds...");
+  //error->warning(FLERR,"Restart functionality not yet tested for granular bonds...");
 }