Skip to content
Browse files

Improved.

  • Loading branch information...
1 parent 0808894 commit ee82d747fe7f0b7b8f8aaacef5637406c8fd08ee @vitroid committed Jun 20, 2012
Showing with 41 additions and 24 deletions.
  1. +3 −3 .depend
  2. +1 −1 .gitignore
  3. +6 −3 Makefile
  4. +3 −1 cluster.f90 → main.f90
  5. +4 −0 physconst.f90
  6. +16 −6 property.f90
  7. +6 −10 settings.f90
  8. +2 −0 test.input
View
6 .depend
@@ -1,7 +1,7 @@
# DO NOT DELETE
-settings.o: property.o vector.o
+main.o: proceed.o property.o force.o settings.o
+settings.o: vector.o
proceed.o: property.o settings.o
-property.o: settings.o vector.o
-cluster.o: proceed.o property.o force.o settings.o
+property.o: physconst.o vector.o
force.o: property.o settings.o
View
2 .gitignore
@@ -1,7 +1,7 @@
@
-cluster
+main
*.o
*.mod
*~
View
9 Makefile
@@ -1,12 +1,15 @@
-SOURCES = vector.f90 proceed.f90 settings.f90 force.f90 cluster.f90 property.f90
+SOURCES = vector.f90 proceed.f90 settings.f90 force.f90 main.f90 property.f90 physconst.f90
OBJECTS = ${patsubst %.f90, %.o, $(SOURCES)}
+TARGET = main
FC = gfortran
-all: cluster
+all: $(TARGET)
%.o: %.f90
$(FC) -c $(FFLAGS) $< -o $@
-cluster: $(OBJECTS)
+main: $(OBJECTS)
$(FC) $(OBJECTS) -o $@
depend:
perl ./f90makedepend $(SOURCES) > .depend
+clean:
+ rm $(TARGET) *~ *.o *.mod
include .depend
View
4 cluster.f90 → main.f90
@@ -14,8 +14,10 @@ program clustermd
call force_to_accel()
call proceed_velocity(dt)
call property_kineticenergy()
- write(STDOUT,*) i,ep,ek,ep+ek
call proceed_position(dt/2)
+ if (logging_interval > 0 .and. mod(i,logging_interval) == 0) then
+ write(STDOUT,*) i,ep,ek,ep+ek, temperature
+ endif
enddo
call settings_write
end program clustermd
View
4 physconst.f90
@@ -0,0 +1,4 @@
+module physconst
+ implicit none
+ real(kind=8), parameter :: kB = 8.314d0
+end module physconst
View
22 property.f90
@@ -1,38 +1,48 @@
module property
use vector
+ use physconst
implicit none
real(kind=8) :: ep,ek ! kJ/mol
+ real(kind=8) :: temperature ! K
+ integer :: num_molecule
type(vector3), allocatable :: position(:) ! Angstrom
type(vector3), allocatable :: velocity(:) ! Angstrom/ps
type(vector3), allocatable :: force_acc(:) ! Force or Acceleration
+ real(kind=8) :: mass ! atomic unit
+ real(kind=8) :: argon_sigma ! Angstrom
+ real(kind=8) :: argon_eps4 ! kJ/mol x 4
contains
- subroutine property_prepare(num_molecule)
- integer, intent(IN) :: num_molecule
+ subroutine property_prepare
allocate(position(num_molecule))
allocate(velocity(num_molecule))
allocate(force_acc(num_molecule))
+ mass = 39.95 !atomic unit
+ argon_sigma = 3.41d0 !Angstrom
+ argon_eps4 = 120d0 * 1d-3 * kB * 4d0 !kJ/mol x 4
end subroutine property_prepare
- subroutine property_read_ar3a(FILE, num_molecule)
- integer, intent(IN) :: FILE, num_molecule
+ subroutine property_read_ar3a(FILE)
+ integer, intent(IN) :: FILE
!local variables
integer :: i
+ read(FILE,*) num_molecule
+ call property_prepare
do i=1,num_molecule
call vector3_read(FILE,position(i))
enddo
end subroutine property_read_ar3a
subroutine property_kineticenergy
- use settings
!local variables
integer :: i
ek = 0d0
do i=1,num_molecule
ek = ek + 0.5 * mass * vector3_inner_product(velocity(i), velocity(i))
enddo
! ek in g/mol * (A/ps)^2 == 10 J/mol
- ek = ek * 0.01 ! kJ/mol
+ ek = ek * 0.01 ! kJ/mol 3NkbT = 2Ek
+ temperature = ek * 1d-3 * 2d0 /(3d0 * num_molecule * kB)
end subroutine property_kineticenergy
end module property
View
16 settings.f90
@@ -4,11 +4,8 @@ module settings
integer, parameter :: STDIN = 5
integer, parameter :: STDOUT = 6
integer :: num_loop
- integer :: num_molecule
+ integer :: logging_interval
real(kind=8) :: dt ! pico seconds
- real(kind=8) :: mass ! atomic unit
- real(kind=8) :: argon_sigma ! Angstrom
- real(kind=8) :: argon_eps4 ! kJ/mol x 4
contains
@@ -18,6 +15,7 @@ subroutine settings_read(FILE)
!local variables
character(len=5) :: tag
integer :: i
+ logging_interval = 10
do
read(FILE,*, END=99) tag
if(tag == "@DTPS")then
@@ -26,16 +24,14 @@ subroutine settings_read(FILE)
if(tag == "@MDLP")then
read(FILE,*) num_loop
endif
+ if(tag == "@NLOG")then
+ read(FILE,*) logging_interval
+ endif
if(tag == "@AR3A")then
- read(FILE,*) num_molecule
- call property_prepare(num_molecule)
- call property_read_ar3a(FILE, num_molecule)
+ call property_read_ar3a(FILE)
endif
enddo
99 continue
- mass = 39.95 !atomic unit
- argon_sigma = 3.41d0 !Angstrom
- argon_eps4 = 120d0 * 0.008314 * 4d0 !kJ/mol x 4
end subroutine settings_read
subroutine settings_write
View
2 test.input
@@ -1,3 +1,5 @@
+@NLOG
+1000
@MDLP
1000000
@DTPS

0 comments on commit ee82d74

Please sign in to comment.
Something went wrong with that request. Please try again.