Subroutine ReGrow( Nlj, Nion, MaxMol, MaxSp, Nmol, &
				   Xlj, Ylj, Zlj, TYPElj, Xion, Yion, Zion, TYPEion, &
				   NSTEPS, MaxSteps, STEPSTART, STEPLENGTH, NTRIALS, &
				   MaxBeads, METHOD, MaxInt, MaxReal, INTPARAM, &
				   REALPARAM, BEADTYPE, BONDSAPART, BoxSize, SPECIES, &
				   LENGTHlj, LENGTHion, STARTlj, STARTion, PROB_SP, &
				   NTemp, BETA, Nham, LNWSTATE, XLRCORR, ELRCORR, &
				   Nljgrs, EPS, SIG, CP, ALP, RMAX, Niongrs, CHARGE, &
				   Alpha, Kmax, Nkvec, KX, KY, KZ, CONST, &
				   EXPX, EXPY, EXPZ, SUMQEXPV, SUMQX, SUMQY, SUMQZ, &
				   ULJBOX, ULJLR, UREAL, UFOURIER, USURF, USELF_CH, &
				   USELF_MOL, ULJ_MOL, UION_MOL, UINTRA, SpID, &
				   Success, f14_lj, f14_ion, Seed )

implicit none

! This subroutine regrows a flexible molecule from StartGrowthStep to the 
! last step and then accepts or rejects the move.

! Nlj is the number of LJ beads in the simulation boxes.
! Nion is the number of ionic beads in the simulation boxes.

integer, intent(inout)									:: Nlj
integer, intent(inout)									:: Nion

! MaxMol is the maximum number of molecules per box.
! MaxSp is the maximum number of Species in the simulation.

integer, intent(in)										:: MaxMol
integer, intent(in)										:: MaxSp

! Nmol is the number of molecules in the simulation boxes.
! Nmol(0, iphase) ==> total number of molecules in that phase.
! Nmol(1:MaxSp, iphase) ==> number of molecules of that species in that phase.

integer, dimension(0:MaxSp), intent(inout)				:: Nmol

! Xlj, Ylj, and Zlj are the coordinates of the LJ beads.
! Xion, Yion, and Zion are the coordinates of the ionic beads.

real, dimension(Nlj), intent(inout)						:: Xlj, Ylj, Zlj
real, dimension(Nion), intent(inout)					:: Xion, Yion, Zion

! TYPElj contains the group identity of each LJ bead.
! TYPEion contains the group identity of each ionic bead.

integer, dimension(Nlj), intent(inout)					:: TYPElj
integer, dimension(Nion), intent(inout)					:: TYPEion

! NSTEPS is the number of configurational bias steps it takes to grow a molecule.
! MaxSteps is the maximum number of CB steps.

integer, dimension(MaxSp), intent(in)					:: NSTEPS
integer, intent(in)										:: MaxSteps

! STEPSTART is the start bead of each CB step.
! STEPLENGTH is the number of beads in each CB step.
! NTRIALS is the number of trial orientations/locations for each CB step.

integer, dimension(MaxSteps, MaxSp), intent(in)			:: STEPSTART
integer, dimension(MaxSteps, MaxSp), intent(in)			:: STEPLENGTH
integer, dimension(MaxSteps, MaxSp), intent(in)			:: NTRIALS

! MaxBeads in the maximum number of beads per molecule.
! METHOD is the method used to grow each bead.
! MaxInt is the maximum number of integer parameters.
! MaxReal is the maximum number of real	parameters.
! INTPARAM are the integer parameters needed by the method to grow each bead.
! REALPARAM are the real parameters needed by the method to grow each bead.
! BEADTYPE indicates whether the bead is 'LJ' or 'ION'.
! BONDSAPART gives the bondlengths any two LJ beads are away from each other.

integer, intent(in)										:: MaxBeads
character*10, dimension(MaxBeads, MaxSp), intent(in)	:: METHOD
integer, intent(in)										:: MaxInt
integer, intent(in)										:: MaxReal
integer, dimension(MaxInt, MaxBeads, MaxSp), intent(in)	:: INTPARAM
real, dimension(MaxReal, MaxBeads, MaxSp), intent(in)	:: REALPARAM
character*5, dimension(MaxBeads, MaxSp), intent(in)		:: BEADTYPE
integer, dimension(MaxBeads,MaxBeads,MaxSp), intent(in)	:: BONDSAPART

! BoxSize is the length of the simulation box.

real, intent(in)										:: BoxSize

! SPECIES contains the species identity of each molecule.
! LENGTHlj contains the number of LJ beads in each molecule.
! LENGTHion	contains the number of ionic beads in each molecule.
! STARTlj contains the starting LJ bead number for each molecule.
! STARTion contains the starting ionic bead number for each molecule.

integer, dimension(Nmol(0)), intent(inout)				:: SPECIES
integer, dimension(Nmol(0)), intent(inout)				:: LENGTHlj, LENGTHion
integer, dimension(Nmol(0)), intent(inout)				:: STARTlj, STARTion

! PROB_SP contains the accumulative probability of selecting a given species.

real, dimension(MaxSp), intent(in)						:: PROB_SP
! NTemp is the number of temperatures being used.

integer, intent(in)										:: NTemp

! BETA contains the reciprical temperature.

real, dimension(NTemp,1), intent(in)					:: BETA

! Nham is the number of hamiltonians being used.

integer, intent(in)										:: Nham

! LNWSTATE contains the weight of each temperature and hamiltonian.

real, dimension(NTemp, Nham), intent(inout)				:: LNWSTATE

! XLRCORR and ELRCORR contain parameters for the long range LJ energy.

real, dimension(605), intent(in)						:: XLRCORR, ELRCORR

! Nljgrs is the number of LJ groups in the system.
! EPS is a rank 3 array containing the eps_ij parameters for each hamiltonian.
! SIG is a rank 3 array containing the sigma_ij parameters for each hamiltonian.
! CP contains the C_ij parameters for each hamiltonian.
! ALP contains the alpha_ij parameters for each hamiltonian.
! RMAX contains the Rmax_ij parameters for each hamiltonian.
									
integer, intent(in)										:: Nljgrs
real, dimension(Nljgrs, Nljgrs, Nham), intent(in)		:: EPS, SIG, CP, ALP, RMAX

! Niongrs is the number of ionic groups in the system.
! CHARGE is a rank 2 array containing the charge of group i for each hamiltonian.
									
integer, intent(in)										:: Niongrs
real, dimension(Niongrs, Nham)							:: CHARGE

! Alpha is an Ewald sum parameter, Alpha = kappa * L, for kappa in A + T.

real, intent(in)										:: Alpha

! Kmax is an Ewald sum parameter.
! Nkvec is the number of k-vectors used in the Fourier sum.
! KX, KY, KZ contain the vector identity of the Nkvec vectors.
! CONST contains the constant part of the Fourier summation for a given Nkvec.

integer, intent(in)										:: Kmax
integer, intent(in)										:: Nkvec
integer, dimension(Nkvec), intent(in)					:: KX, KY, KZ
real, dimension(Nkvec), intent(in)						:: CONST

! EXPX contains the value of exp( i*kx*x ) for a given kx and ion.

complex, dimension(0:Kmax, Nion), intent(inout)			:: EXPX
complex, dimension(-Kmax:Kmax, Nion), intent(inout)		:: EXPY, EXPZ

! SUMQEXPV contains the summation of qi*exp(i*(kx*x + ky*y + kz*z)) 
! for a given k-vector and hamiltonian.

complex, dimension(Nkvec, Nham), intent(inout)			:: SUMQEXPV

! SUMQX is the summation of qi * xi for all ions in the box.

real, dimension(Nham), intent(inout)					:: SUMQX, SUMQY, SUMQZ

! ULJ is the LJ energy of the system without the long range correction.
! UFOURIER is the coulombic fourier energy of the system.
! UREAL is the coulombic real energy of the system.
! USURF is the coulombic surface energy of the system.
! USELF_CH is the summation of the square of all the charges.
! USELF_MOL is the self energy of a given molecule.
! ULJ_MOL is the non bonded intramolecular LJ energy of a molecule.
! UION_MOL is the non bonded intramolecular ionic energy of a molecule.
! UINTRA contains the bending and torsional energy of each molecule.

real, dimension(Nham), intent(inout)					:: ULJBOX, ULJLR
real, dimension(Nham), intent(inout)					:: UREAL, UFOURIER
real, dimension(Nham), intent(inout)					:: USURF, USELF_CH
real, dimension(MaxMol,Nham), intent(inout)				:: USELF_MOL
real, dimension(MaxMol,Nham), intent(inout)				:: ULJ_MOL
real, dimension(MaxMol,Nham), intent(inout)				:: UION_MOL
real, dimension(MaxMol), intent(inout)					:: UINTRA

! SpID indicates which type of species was used in the move.

integer, intent(out)									:: SpID

! Success is a logical indicating whether the move was successful or not.

logical, intent(out)									:: Success

! f14_lj and f14_ion are damping factors for the 1-4 intramolecular interactions

real, intent(in)										:: f14_lj
real, intent(in)										:: f14_ion

! Seed is the current random number generator seed value.

integer, intent(inout)									:: Seed

real, external											:: ran2

! Local Variables

integer										:: i, k
integer										:: Count
integer										:: MolSpecies
integer										:: Mol
integer										:: lenlj, lenion
integer										:: lenljst, lenionst
integer										:: stlj, stion
integer										:: endlj, endion
integer										:: StartGrowthStep
integer, dimension(Nlj)						:: TYPElj_new
integer, dimension(Nion)					:: TYPEion_new
integer, dimension(Nljgrs)					:: NGROUPS

integer, allocatable, dimension(:) 			:: TYPElj_grow
integer, allocatable, dimension(:) 			:: TYPEion_grow

logical										:: New

real										:: Random
real										:: PiRatio
real										:: Uintra_grow
real										:: Uintra_new
real										:: CoulCombo
real										:: BigW_old, BigW_new
real, dimension(Nlj)						:: Xlj_new, Ylj_new, Zlj_new
real, dimension(Nion)						:: Xion_new, Yion_new, Zion_new
real, dimension(Nham)						:: ULJBOX_new, ULJLR_new
real, dimension(Nham)						:: UREAL_new, UFOURIER_new
real, dimension(Nham)						:: USURF_new, USELF_CH_new
real, dimension(Nham)						:: USELF_MOL_new
real, dimension(Nham)						:: UION_MOL_new
real, dimension(Nham)						:: ULJ_MOL_new
real, dimension(Nham)						:: ULJBOX_grow, ULJLR_grow
real, dimension(Nham)						:: UREAL_grow, UFOURIER_grow
real, dimension(Nham)						:: USURF_grow, USELF_CH_grow
real, dimension(Nham)						:: USELF_MOL_grow
real, dimension(Nham)						:: UION_MOL_grow
real, dimension(Nham)						:: ULJ_MOL_grow
real, dimension(Nham)						:: dUFOURIER
real, dimension(Nham)						:: dUSURF
real, dimension(Nham)						:: SUMQX_new, SUMQY_new, SUMQZ_new
real, dimension(Nham)						:: SUMQX_grow, SUMQY_grow, SUMQZ_grow
real, dimension(Nljgrs,Nljgrs,Nham)			:: EPS_TMP

real, allocatable, dimension(:) 			:: Xlj_grow1, Ylj_grow1, Zlj_grow1
real, allocatable, dimension(:) 			:: Xion_grow1, Yion_grow1, Zion_grow1
real, allocatable, dimension(:) 			:: Xlj_grow2, Ylj_grow2, Zlj_grow2
real, allocatable, dimension(:) 			:: Xion_grow2, Yion_grow2, Zion_grow2

real, parameter								:: Pi = 3.14159265359
real, parameter								:: ec = 1.60217733e-19
real, parameter								:: eps0 = 8.854187817e-12
real, parameter								:: kB = 1.380658e-23
										
complex, dimension(Nkvec, Nham)				:: SUMQEXPV_new, SUMQEXPV_grow
complex, allocatable, dimension(:,:)		:: EXPX_grow1, EXPY_grow1, EXPZ_grow1
complex, allocatable, dimension(:,:)		:: EXPX_grow2, EXPY_grow2, EXPZ_grow2
complex, allocatable, dimension(:,:)		:: CZERO1, CZERO2


CoulCombo = ec * ec * 1.0e10 / ( 4.0 * Pi * eps0 * kB )

Success = .False.

Random = ran2(Seed)
SpID = 0
i = 1

do while ( SpID == 0 )

	if( Random < PROB_SP(i) ) SpID = i
	
	i = i + 1

end do

if( Nmol(SpID) == 0 ) return

MolSpecies = int( Nmol( SpID) * ran2(Seed) ) + 1

i = 0
Count = 0

do while ( Count < MolSpecies )
	
	i = i + 1
	
	if( SPECIES(i) == SpID ) Count = Count + 1

end do

Mol = i

lenlj = LENGTHlj(Mol)
stlj = STARTlj(Mol)
endlj = stlj + lenlj - 1

lenion = LENGTHion(Mol)
stion = STARTion(Mol)
endion = stion + lenion - 1

if( stlj == 1 ) then

	if( Nmol(0) == 1 ) then

		Xlj_new = 0.0
		Ylj_new = 0.0
		Zlj_new = 0.0

		TYPElj_new = 0
	
	else

		Xlj_new( 1:Nlj-lenlj ) = Xlj( endlj+1:Nlj )
		Ylj_new( 1:Nlj-lenlj ) = Ylj( endlj+1:Nlj )
		Zlj_new( 1:Nlj-lenlj ) = Zlj( endlj+1:Nlj )

		TYPElj_new( 1:Nlj-lenlj ) = TYPElj( endlj+1:Nlj )

	end if

else if( stlj + lenlj - 1 == Nlj ) then

	Xlj_new( 1:stlj-1 ) = Xlj( 1:stlj-1 )
	Ylj_new( 1:stlj-1 ) = Ylj( 1:stlj-1 )
	Zlj_new( 1:stlj-1 ) = Zlj( 1:stlj-1 )

	TYPElj_new( 1:stlj-1 ) = TYPElj( 1:stlj-1 )

else

	Xlj_new( 1:stlj-1 ) = Xlj( 1:stlj-1 )
	Xlj_new( stlj:Nlj-lenlj ) = Xlj( endlj+1:Nlj )
	Ylj_new( 1:stlj-1 ) = Ylj( 1:stlj-1 )
	Ylj_new( stlj:Nlj-lenlj ) = Ylj( endlj+1:Nlj )
	Zlj_new( 1:stlj-1 ) = Zlj( 1:stlj-1 )
	Zlj_new( stlj:Nlj-lenlj ) = Zlj( endlj+1:Nlj )

	TYPElj_new( 1:stlj-1 ) = TYPElj( 1:stlj-1 )
	TYPElj_new( stlj:Nlj-lenlj ) = TYPElj( endlj+1:Nlj )

end if

if( lenion > 0 ) then

	if( stion == 1 ) then

		if( Nmol(0) == 1 ) then

			Xion_new = 0.0
			Yion_new = 0.0
			Zion_new = 0.0

			TYPEion_new = 0

		else

			Xion_new( 1:Nion-lenion ) = Xion( endion+1:Nion )
			Yion_new( 1:Nion-lenion ) = Yion( endion+1:Nion )
			Zion_new( 1:Nion-lenion ) = Zion( endion+1:Nion )

			TYPEion_new( 1:Nion-lenion ) = TYPEion( endion+1:Nion )

		end if

	else if( stion + lenion - 1 == Nion ) then

		Xion_new( 1:stion-1 ) = Xion( 1:stion-1 )
		Yion_new( 1:stion-1 ) = Yion( 1:stion-1 )
		Zion_new( 1:stion-1 ) = Zion( 1:stion-1 )

		TYPEion_new( 1:stion-1 ) = TYPEion( 1:stion-1 )

	else

		Xion_new( 1:stion-1 ) = Xion( 1:stion-1 )
		Xion_new( stion:Nion-lenion ) = Xion( endion+1:Nion )
		Yion_new( 1:stion-1 ) = Yion( 1:stion-1 )
		Yion_new( stion:Nion-lenion ) = Yion( endion+1:Nion )
		Zion_new( 1:stion-1 ) = Zion( 1:stion-1 )
		Zion_new( stion:Nion-lenion ) = Zion( endion+1:Nion )

		TYPEion_new( 1:stion-1 ) = TYPEion( 1:stion-1 )
		TYPEion_new( stion:Nion-lenion ) = TYPEion( endion+1:Nion )

	end if

end if

allocate( Xlj_grow1(lenlj) )
allocate( Ylj_grow1(lenlj) )
allocate( Zlj_grow1(lenlj) )

allocate( Xlj_grow2(lenlj) )
allocate( Ylj_grow2(lenlj) )
allocate( Zlj_grow2(lenlj) )

allocate( TYPElj_grow(lenlj) )

! For n-alkanes only.

!if( ran2(seed) < 0.5 ) then

!	do i = 1, lenlj
 
!		Xlj_grow1(i) = Xlj( endlj - i + 1 )
!		Ylj_grow1(i) = Ylj( endlj - i + 1 )
!		Zlj_grow1(i) = Zlj( endlj - i + 1 )
 
!	end do 
 
!else 
 
	Xlj_grow1( 1:lenlj ) = Xlj( stlj:endlj )
	Ylj_grow1( 1:lenlj ) = Ylj( stlj:endlj )
	Zlj_grow1( 1:lenlj ) = Zlj( stlj:endlj )
 
!end if

Xlj_grow2 = Xlj_grow1
Ylj_grow2 = Ylj_grow1
Zlj_grow2 = Zlj_grow1

TYPElj_grow( 1:lenlj ) = TYPElj( stlj:endlj )

!StartGrowthStep = 1

StartGrowthStep = int( ran2(Seed) * ( NSTEPS(SpID) - 1 ) ) + 2

! For n-alkanes only

!StartGrowthStep = NSTEPS(SpID) / 2 + int( ran2(Seed) * ( NSTEPS(SpID) / 2 ) ) + 1 

lenljst = 0
lenionst = 0

do i = StartGrowthStep, NSTEPS(SpID)

	do k = STEPSTART(i,SpID), STEPSTART(i,SpID) + STEPLENGTH(i,SpID) - 1
	
		if(	BEADTYPE(k,SpID) == 'LJ' ) then

			lenljst = lenljst + 1
		
		else if( BEADTYPE(k,SpID) == 'ION' ) then

			lenionst = lenionst + 1

		end if

	end do

end do

ULJLR_new = ULJLR
ULJLR_grow = ULJLR

if( lenljst > 0 ) then

	NGROUPS = 0
		
	do k = 1, Nlj-lenlj

		NGROUPS( TYPElj_new(k) ) = NGROUPS( TYPElj_new(k) ) + 1

	end do

	do k = 1, lenlj - lenljst

		NGROUPS( TYPElj_grow(k) ) = NGROUPS( TYPElj_grow(k) ) + 1

	end do

	if( CP(1,1,1) > 0.0 ) then

		EPS_TMP = EPS * CP * ALP / ( ALP - 6.0 ) / 4.0

	else

		EPS_TMP = EPS

	end if

	call lrcorr( Nljgrs, Nham, BoxSize, NGROUPS, EPS_TMP, SIG, XLRCORR, &
				 ELRCORR, ULJLR_new )

	ULJLR_grow = ULJLR_new

end if

ULJBOX_new = 0.0
ULJBOX_grow = 0.0

ULJ_MOL_new = 0.0
ULJ_MOL_grow = 0.0

UION_MOL_new = 0.0
UION_MOL_grow = 0.0

Uintra_new = 0.0
Uintra_grow = 0.0

if( lenionst > 0 ) then

	allocate( Xion_grow1(lenion) )
	allocate( Yion_grow1(lenion) )
	allocate( Zion_grow1(lenion) )

	allocate( Xion_grow2(lenion) )
	allocate( Yion_grow2(lenion) )
	allocate( Zion_grow2(lenion) )

	allocate( TYPEion_grow(lenion) )

	Xion_grow1( 1:lenion ) = Xion( stion:endion )
	Yion_grow1( 1:lenion ) = Yion( stion:endion )
	Zion_grow1( 1:lenion ) = Zion( stion:endion )

	Xion_grow2 = Xion_grow1
	Yion_grow2 = Yion_grow1
	Zion_grow2 = Zion_grow1

	TYPEion_grow( 1:lenion ) = TYPEion( stion:endion )

	UREAL_new = 0.0
	UREAL_grow = 0.0

	Xion_grow1 = -Xion_grow1
	Yion_grow1 = -Yion_grow1
	Zion_grow1 = -Zion_grow1

	call Surf_Move( lenionst, Xion_grow1(lenion-lenionst+1:lenion), &
					Yion_grow1(lenion-lenionst+1:lenion), &
					Zion_grow1(lenion-lenionst+1:lenion), &
					TYPEion_grow(lenion-lenionst+1:lenion), &
					Nham, Niongrs, CHARGE, &
					BoxSize, SUMQX, SUMQY, SUMQZ, &
					SUMQX_new, SUMQY_new, SUMQZ_new, dUSURF )

	Xion_grow1 = -Xion_grow1
	Yion_grow1 = -Yion_grow1
	Zion_grow1 = -Zion_grow1

	SUMQX_grow = SUMQX_new
	SUMQY_grow = SUMQY_new
	SUMQZ_grow = SUMQZ_new


	USURF_new = USURF + dUSURF * CoulCombo

	USURF_grow = USURF_new

	allocate( EXPX_grow1(0:Kmax, lenion ) )
	allocate( EXPY_grow1(-Kmax:Kmax, lenion ) )
	allocate( EXPZ_grow1(-Kmax:Kmax, lenion ) )

	allocate( EXPX_grow2(0:Kmax, lenion ) )
	allocate( EXPY_grow2(-Kmax:Kmax, lenion ) )
	allocate( EXPZ_grow2(-Kmax:Kmax, lenion ) )

	allocate( CZERO1(0:Kmax,lenion) )
	allocate( CZERO2(-Kmax:Kmax,lenion) )

	CZERO1 = ( 0.0, 0.0 )
	CZERO2 = ( 0.0, 0.0 )

	EXPX_grow1(:, 1:lenion ) = EXPX( :, stion:endion )
	EXPY_grow1(:, 1:lenion ) = EXPY( :, stion:endion )
	EXPZ_grow1(:, 1:lenion ) = EXPZ( :, stion:endion )

	EXPX_grow2 = EXPX_grow1
	EXPY_grow2 = EXPY_grow1
	EXPZ_grow2 = EXPZ_grow1

	CHARGE = -CHARGE

	call Fourier_Move( lenionst, Xion_grow1(lenion-lenionst+1:lenion), &
					   Yion_grow1(lenion-lenionst+1:lenion), &
					   Zion_grow1(lenion-lenionst+1:lenion), &
					   TYPEion_grow(lenion-lenionst+1:lenion), &
					   Nham, Niongrs, CHARGE, BoxSize, &
					   Kmax, Nkvec, KX, KY, KZ, &
					   CONST, CZERO1, CZERO2, CZERO2, &
					   EXPX_grow1(:,lenion-lenionst+1:lenion), &
					   EXPY_grow1(:,lenion-lenionst+1:lenion), &
					   EXPZ_grow1(:,lenion-lenionst+1:lenion), &
					   SUMQEXPV, SUMQEXPV_new, dUFOURIER )

	CHARGE = -CHARGE

	SUMQEXPV_grow = SUMQEXPV_new

	UFOURIER_new = UFOURIER + dUFOURIER * CoulCombo

	UFOURIER_grow = UFOURIER_new

	deallocate( CZERO1 )
	deallocate( CZERO2 )

	USELF_CH_new = 0.0
	USELF_CH_grow = 0.0

	call SelfMolecule( lenion-lenionst, Xion_grow1, Yion_grow1, Zion_grow1, &
					   TYPEion_grow, Nham, Niongrs, CHARGE, BoxSize, &
					   Alpha, USELF_MOL_new )

	USELF_MOL_new = USELF_MOL_new * CoulCombo
	USELF_MOL_grow = USELF_MOL_new

else

	UREAL_new = 0.0
	UREAL_grow = 0.0
	USURF_new = 0.0
	USURF_grow = 0.0
	UFOURIER_new = 0.0
	UFOURIER_grow = 0.0
	USELF_CH_new = 0.0
	USELF_CH_grow = 0.0
	USELF_MOL_new = 0.0
	USELF_MOL_grow = 0.0
	
end if

New = .False.

call Grow( NSTEPS(SpID), STEPSTART(:,SpID), STEPLENGTH(:,SpID), &
		   NTRIALS(:,SpID), lenlj, lenion, MaxBeads, METHOD(:,SpID), &
		   MaxInt, MaxReal, INTPARAM(:,:,SpID), REALPARAM(:,:,SpID), &
		   BEADTYPE(:,SpID), New, StartGrowthStep, Xlj_grow1, Ylj_grow1, &
		   Zlj_grow1, TYPElj_grow, Xion_grow1, Yion_grow1, &
		   Zion_grow1, TYPEion_grow, Nham, ULJBOX_new, &
		   ULJLR_new, ULJ_MOL_new, UREAL_new, USURF_new, &
		   UFOURIER_new, USELF_CH_new, USELF_MOL_new, UION_MOL_new, &
		   Uintra_new, Niongrs, CHARGE, Alpha, Kmax, Nkvec, KX, KY, KZ, &
		   CONST, SUMQX_new, SUMQY_new, SUMQZ_new, EXPX_grow1, EXPY_grow1, &
		   EXPZ_grow1, SUMQEXPV_new, BigW_old, NTemp, BETA, LNWSTATE, &
		   Nlj-lenlj, Xlj_new, Ylj_new, Zlj_new, TYPElj_new, &
		   Nion-lenion, Xion_new, Yion_new, Zion_new, TYPEion_new, &
		   BoxSize, Nljgrs, EPS, SIG, CP, ALP, RMAX, XLRCORR, ELRCORR, &
		   BONDSAPART(:,:,SpID), f14_lj, f14_ion, Seed )


New = .True.

call Grow( NSTEPS(SpID), STEPSTART(:,SpID), STEPLENGTH(:,SpID), &
		   NTRIALS(:,SpID), lenlj, lenion, MaxBeads, METHOD(:,SpID), &
		   MaxInt, MaxReal, INTPARAM(:,:,SpID), REALPARAM(:,:,SpID), &
		   BEADTYPE(:,SpID), New, StartGrowthStep, Xlj_grow2, Ylj_grow2, &
		   Zlj_grow2, TYPElj_grow, Xion_grow2, Yion_grow2, &
		   Zion_grow2, TYPEion_grow, Nham, ULJBOX_grow, &
		   ULJLR_grow, ULJ_MOL_grow, UREAL_grow, USURF_grow, &
		   UFOURIER_grow, USELF_CH_grow, USELF_MOL_grow, UION_MOL_grow, &
		   Uintra_grow, Niongrs, CHARGE, Alpha, Kmax, Nkvec, KX, KY, KZ, &
		   CONST, SUMQX_grow, SUMQY_grow, SUMQZ_grow, EXPX_grow2, EXPY_grow2, &
		   EXPZ_grow2, SUMQEXPV_grow, BigW_new, NTemp, BETA, LNWSTATE, &
		   Nlj-lenlj, Xlj_new, Ylj_new, Zlj_new, TYPElj_new, &
		   Nion-lenion, Xion_new, Yion_new, Zion_new, TYPEion_new, &
		   BoxSize, Nljgrs, EPS, SIG, CP, ALP, RMAX, XLRCORR, ELRCORR, &
		   BONDSAPART(:,:,SpID), f14_lj, f14_ion, Seed )

if( BigW_new < 1.0e-250 ) then

	PiRatio = -2

else

	PiRatio = BigW_new / BigW_old

end if

if( ran2(Seed) < PiRatio ) then

	Success = .True.
	
	ULJBOX_grow = ULJBOX + ULJBOX_grow - ULJBOX_new
	ULJ_MOL_grow(:) = ULJ_MOL(Mol,:) + ULJ_MOL_grow(:) - ULJ_MOL_new(:)

	UREAL_grow = UREAL + UREAL_grow - UREAL_new
	USELF_CH_grow = USELF_CH
	UION_MOL_grow(:) = UION_MOL(Mol,:) + UION_MOL_grow(:) - UION_MOL_new(:)

	Uintra_grow	= UINTRA(Mol) + Uintra_grow - Uintra_new

	Xlj( stlj:endlj ) = Xlj_grow2( 1:lenlj )
	Ylj( stlj:endlj ) = Ylj_grow2( 1:lenlj )
	Zlj( stlj:endlj ) = Zlj_grow2( 1:lenlj )

	ULJBOX = ULJBOX_grow
	ULJLR = ULJLR_grow

	if( lenion > 0 ) then

		Xion( stion:endion ) = Xion_grow2( 1:lenion )
		Yion( stion:endion ) = Yion_grow2( 1:lenion )
		Zion( stion:endion ) = Zion_grow2( 1:lenion )

		UREAL = UREAL_grow
		USURF = USURF_grow
		UFOURIER = UFOURIER_grow

		SUMQX = SUMQX_grow
		SUMQY = SUMQY_grow
		SUMQZ = SUMQZ_grow

		SUMQEXPV = SUMQEXPV_grow

		EXPX( :, stion:endion ) = EXPX_grow2( :, 1:lenion )
		EXPY( :, stion:endion ) = EXPY_grow2( :, 1:lenion )
		EXPZ( :, stion:endion ) = EXPZ_grow2( :, 1:lenion )

	end if

	USELF_MOL(Mol,:) = USELF_MOL_grow(:)

	UION_MOL(Mol,:) = UION_MOL_grow(:)

	ULJ_MOL(Mol,:) = ULJ_MOL_grow(:)

	UINTRA(Mol) = Uintra_grow

end if

deallocate( Xlj_grow1 )
deallocate( Ylj_grow1 )
deallocate( Zlj_grow1 )

deallocate( Xlj_grow2 )
deallocate( Ylj_grow2 )
deallocate( Zlj_grow2 )

deallocate( TYPElj_grow )

if( lenion > 0 ) then

	deallocate( Xion_grow1 )
	deallocate( Yion_grow1 )
	deallocate( Zion_grow1 )

	deallocate( Xion_grow2 )
	deallocate( Yion_grow2 )
	deallocate( Zion_grow2 )

	deallocate( TYPEion_grow )

	deallocate( EXPX_grow1 )
	deallocate( EXPY_grow1 )
	deallocate( EXPZ_grow1 )

	deallocate( EXPX_grow2 )
	deallocate( EXPY_grow2 )
	deallocate( EXPZ_grow2 )

end if

return

end subroutine ReGrow