dcdfort
Modern Fortran library for analyzing DCD trajectories
Data Types | Functions/Subroutines | Variables
dcdfort_writer Module Reference

Data Types

type  dcdwriter
 

Functions/Subroutines

subroutine dcdwriter_open (this, filename)
 
subroutine dcdwriter_write_header (this, istart, nevery, timestep, natoms)
 
subroutine dcdwriter_close (this)
 
subroutine dcdwriter_write_next (this, xyz, box_in)
 

Variables

real(8), parameter pi = 2.0d0*dacos(0.0d0)
 

Function/Subroutine Documentation

◆ dcdwriter_close()

subroutine dcdfort_writer::dcdwriter_close ( class(dcdwriter), intent(inout)  this)
124 
125  implicit none
126  class(dcdwriter), intent(inout) :: this
127 
128  close(this%u)
129 

◆ dcdwriter_open()

subroutine dcdfort_writer::dcdwriter_open ( class(dcdwriter), intent(inout)  this,
character (len=*)  filename 
)
43 
44  implicit none
45  character (len=*) :: filename
46  class(dcdwriter), intent(inout) :: this
47  integer :: filesize
48  logical :: ex
49 
50  open(newunit=this%u, file=trim(filename), form="unformatted", access="stream", status="new")
51 

◆ dcdwriter_write_header()

subroutine dcdfort_writer::dcdwriter_write_header ( class(dcdwriter), intent(inout)  this,
integer  istart,
integer  nevery,
real  timestep,
integer  natoms 
)
55 
56  implicit none
57 
58  integer :: dummy, istart, nevery, natoms, i, n
59  real :: timestep
60  class(dcdwriter), intent(inout) :: this
61  character (len=79) :: remarks1, remarks2
62  character (len=8) :: date
63  character (len=10) :: time
64 
65  call date_and_time(date=date,time=time)
66  remarks1 = "Created by libdcdfort"
67  remarks2 = "REMARK Created on "//date//" "//time
68 
69  this%nframes = 0
70  this%iend = istart
71  this%nevery = nevery
72 
73  write(this%u) 84
74  write(this%u) "CORD"
75 
76  inquire(unit=this%u, pos=this%nframes_pos)
77  ! Number of snapshots in file
78  write(this%u) this%nframes
79 
80  ! Timestep of first snapshot
81  write(this%u) istart
82 
83  ! Save snapshots every this many steps
84  write(this%u) nevery
85 
86  inquire(unit=this%u, pos=this%iend_pos)
87  ! Timestep of last snapshot
88  write(this%u) this%iend
89 
90  do i = 1, 4
91  write(this%u) 0
92  end do
93 
94  ! Has unit cell
95  write(this%u) 1
96 
97  ! Simulation timestep
98  write(this%u) timestep
99 
100  do i = 1, 9
101  write(this%u) 0
102  end do
103 
104  ! Pretend to be CHARMM version 24
105  write(this%u) 24
106  write(this%u) 84
107  write(this%u) 164
108 
109  write(this%u) 2
110  write(this%u) remarks1//c_null_char
111  write(this%u) remarks2//c_null_char
112 
113  write(this%u) 164
114  write(this%u) 4
115 
116  ! Number of atoms in each snapshot
117  write(this%u) natoms
118 
119  write(this%u) 4
120 

◆ dcdwriter_write_next()

subroutine dcdfort_writer::dcdwriter_write_next ( class(dcdwriter), intent(inout)  this,
real, dimension(:,:), intent(in)  xyz,
real(8), dimension(6), intent(in)  box_in 
)
133 
134  implicit none
135  real, intent(in) :: xyz(:,:)
136  real(8), intent(in) :: box_in(6)
137  real(8) :: box(6)
138  integer :: dummy
139  class(dcdwriter), intent(inout) :: this
140  integer :: coord_size
141 
142  coord_size = size(xyz,2)*4
143  box = box_in
144 
145  ! Should be 48 (6 double precision floats)
146  write(this%u) 48
147 
148  box(4) = (90.0 - box(4)) * dsin(pi/180.0)
149  box(5) = (90.0 - box(5)) * dsin(pi/180.0)
150  box(6) = (90.0 - box(6)) * dsin(pi/180.0)
151  write(this%u) box(1) ! A
152  write(this%u) box(6) ! gamma
153  write(this%u) box(2) ! B
154  write(this%u) box(5) ! beta
155  write(this%u) box(4) ! alpha
156  write(this%u) box(3) ! C
157 
158  write(this%u) 48
159  write(this%u) coord_size
160 
161  write(this%u) xyz(1,:)
162 
163  write(this%u) 48
164  write(this%u) coord_size
165 
166  write(this%u) xyz(2,:)
167 
168  write(this%u) 48
169  write(this%u) coord_size
170 
171  write(this%u) xyz(3,:)
172 
173  write(this%u) 48
174 
175  inquire(unit=this%u, pos=this%curr_pos)
176 
177  this%nframes = this%nframes+1
178  this%iend = this%iend + this%nevery
179 
180  write(this%u, pos=this%nframes_pos) this%nframes
181  write(this%u, pos=this%iend_pos) this%iend
182  write(this%u, pos=this%curr_pos)
183 

Variable Documentation

◆ pi

real(8), parameter dcdfort_writer::pi = 2.0d0*dacos(0.0d0)
28  real(8), parameter :: pi = 2.0d0*dacos(0.0d0)