hdiff output

r32617/cost_function.cu 2017-05-26 12:30:26.352901143 +0100 r32616/cost_function.cu 2017-05-26 12:30:28.728932123 +0100
 87:  87: 
 88:         CudaSafeCall( cudaMalloc(&m_d_xRigid,  m_nDegFreedom * sizeof(double)) ); 88:         CudaSafeCall( cudaMalloc(&m_d_xRigid,  m_nDegFreedom * sizeof(double)) );
 89:         CudaSafeCall( cudaMalloc(&m_d_gkRigid, m_nDegFreedom * sizeof(double)) ); 89:         CudaSafeCall( cudaMalloc(&m_d_gkRigid, m_nDegFreedom * sizeof(double)) );
 90:  90: 
 91:         if (m_coords != NULL) { 91:         if (m_coords != NULL) {
 92:                 CudaSafeCall( cudaMalloc(&m_d_coords,   numDimensions * sizeof(double)) ); 92:                 CudaSafeCall( cudaMalloc(&m_d_coords,   numDimensions * sizeof(double)) );
 93:                 CudaSafeCall( cudaMalloc(&m_d_energy, sizeof(double)) ); 93:                 CudaSafeCall( cudaMalloc(&m_d_energy, sizeof(double)) );
 94:                 CudaSafeCall( cudaMalloc(&m_d_gradient, numDimensions * sizeof(double)) ); 94:                 CudaSafeCall( cudaMalloc(&m_d_gradient, numDimensions * sizeof(double)) );
 95:         } 95:         }
 96:  96: 
 97:         if (m_isAtomFrozen != NULL) { 97:         CudaSafeCall( cudaMalloc(&m_d_isAtomFrozen, numDimensions/3 * sizeof(bool)) );
 98:                 CudaSafeCall( cudaMalloc(&m_d_isAtomFrozen, numDimensions/3 * sizeof(bool)) ); 
 99:         } 
100:  98: 
101:         CudaSafeCall( cudaMalloc(&m_d_nLargeRigidBody, sizeof(int)) ); 99:         CudaSafeCall( cudaMalloc(&m_d_nLargeRigidBody, sizeof(int)) );
102:         CudaSafeCall( cudaMalloc(&m_d_largeRigidIndices, m_nLargeRigidBody * sizeof(int)) );100:         CudaSafeCall( cudaMalloc(&m_d_largeRigidIndices, m_nLargeRigidBody * sizeof(int)) );
103: 101: 
104:         CudaSafeCall( cudaMalloc(&m_d_ljAddRep, m_nAddTarget * m_nAddTarget * sizeof(double)) );102:         CudaSafeCall( cudaMalloc(&m_d_ljAddRep, m_nAddTarget * m_nAddTarget * sizeof(double)) );
105:         CudaSafeCall( cudaMalloc(&m_d_ljAddAtt, m_nAddTarget * m_nAddTarget * sizeof(double)) );103:         CudaSafeCall( cudaMalloc(&m_d_ljAddAtt, m_nAddTarget * m_nAddTarget * sizeof(double)) );
106: 104: 
107:         int counter = 0;105:         int counter = 0;
108:         while (counter < m_nLargeRigidBody) {106:         while (counter < m_nLargeRigidBody) {
109:                 for (int i = 0; i < m_nRigidBody; ++i) {107:                 for (int i = 0; i < m_nRigidBody; ++i) {
147:         CudaSafeCall( cudaMemcpy(m_d_rigidSingles, m_rigidSingles, (m_nDegFreedom/3 - 2*m_nRigidBody) * sizeof(int), cudaMemcpyHostToDevice) );145:         CudaSafeCall( cudaMemcpy(m_d_rigidSingles, m_rigidSingles, (m_nDegFreedom/3 - 2*m_nRigidBody) * sizeof(int), cudaMemcpyHostToDevice) );
148:         CudaSafeCall( cudaMemcpy(m_d_rigidInverse, m_rigidInverse, m_nRigidBody * 3 * 3 * sizeof(double), cudaMemcpyHostToDevice) );146:         CudaSafeCall( cudaMemcpy(m_d_rigidInverse, m_rigidInverse, m_nRigidBody * 3 * 3 * sizeof(double), cudaMemcpyHostToDevice) );
149: 147: 
150:         CudaSafeCall( cudaMemcpy(m_d_nRigidBody, &m_nRigidBody, sizeof(int), cudaMemcpyHostToDevice) );148:         CudaSafeCall( cudaMemcpy(m_d_nRigidBody, &m_nRigidBody, sizeof(int), cudaMemcpyHostToDevice) );
151:         CudaSafeCall( cudaMemcpy(m_d_rigidMaxSite, &m_rigidMaxSite, sizeof(int), cudaMemcpyHostToDevice) );149:         CudaSafeCall( cudaMemcpy(m_d_rigidMaxSite, &m_rigidMaxSite, sizeof(int), cudaMemcpyHostToDevice) );
152: 150: 
153:         if (m_coords != NULL) {151:         if (m_coords != NULL) {
154:                 CudaSafeCall( cudaMemcpy(m_d_coords, m_coords, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );152:                 CudaSafeCall( cudaMemcpy(m_d_coords, m_coords, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );
155:         }153:         }
156: 154: 
157:         if (m_isAtomFrozen != NULL) {155:         CudaSafeCall( cudaMemcpy(m_d_isAtomFrozen, m_isAtomFrozen, numDimensions/3 * sizeof(bool), cudaMemcpyHostToDevice) );
158:                 CudaSafeCall( cudaMemcpy(m_d_isAtomFrozen, m_isAtomFrozen, numDimensions/3 * sizeof(bool), cudaMemcpyHostToDevice) ); 
159:         } 
160: 156: 
161:         CudaSafeCall( cudaMemcpy(m_d_zerosx, zerosx, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );157:         CudaSafeCall( cudaMemcpy(m_d_zerosx, zerosx, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );
162:         CudaSafeCall( cudaMemcpy(m_d_zerosy, zerosy, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );158:         CudaSafeCall( cudaMemcpy(m_d_zerosy, zerosy, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );
163:         CudaSafeCall( cudaMemcpy(m_d_zerosz, zerosz, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );159:         CudaSafeCall( cudaMemcpy(m_d_zerosz, zerosz, numDimensions * sizeof(double), cudaMemcpyHostToDevice) );
164: 160: 
165:         CudaSafeCall( cudaMemcpy(m_d_nLargeRigidBody, &m_nLargeRigidBody, sizeof(int), cudaMemcpyHostToDevice) );161:         CudaSafeCall( cudaMemcpy(m_d_nLargeRigidBody, &m_nLargeRigidBody, sizeof(int), cudaMemcpyHostToDevice) );
166:         CudaSafeCall( cudaMemcpy(m_d_largeRigidIndices, m_largeRigidIndices, m_nLargeRigidBody * sizeof(int), cudaMemcpyHostToDevice) );162:         CudaSafeCall( cudaMemcpy(m_d_largeRigidIndices, m_largeRigidIndices, m_nLargeRigidBody * sizeof(int), cudaMemcpyHostToDevice) );
167: 163: 
168:         CudaSafeCall( cudaMemcpy(m_d_ljAddRep, m_ljAddRep, m_nAddTarget * m_nAddTarget * sizeof(double), cudaMemcpyHostToDevice) );164:         CudaSafeCall( cudaMemcpy(m_d_ljAddRep, m_ljAddRep, m_nAddTarget * m_nAddTarget * sizeof(double), cudaMemcpyHostToDevice) );
169:         CudaSafeCall( cudaMemcpy(m_d_ljAddAtt, m_ljAddAtt, m_nAddTarget * m_nAddTarget * sizeof(double), cudaMemcpyHostToDevice) );165:         CudaSafeCall( cudaMemcpy(m_d_ljAddAtt, m_ljAddAtt, m_nAddTarget * m_nAddTarget * sizeof(double), cudaMemcpyHostToDevice) );
196: 192: 
197:         CudaSafeCall( cudaFree(m_d_nRigidBody) );193:         CudaSafeCall( cudaFree(m_d_nRigidBody) );
198:         CudaSafeCall( cudaFree(m_d_rigidMaxSite) );194:         CudaSafeCall( cudaFree(m_d_rigidMaxSite) );
199: 195: 
200:         if (m_coords != NULL) {196:         if (m_coords != NULL) {
201:                 CudaSafeCall( cudaFree(m_d_coords) );197:                 CudaSafeCall( cudaFree(m_d_coords) );
202:                 CudaSafeCall( cudaFree(m_d_energy) );198:                 CudaSafeCall( cudaFree(m_d_energy) );
203:                 CudaSafeCall( cudaFree(m_d_gradient) );199:                 CudaSafeCall( cudaFree(m_d_gradient) );
204:         }200:         }
205: 201: 
206:         if (m_isAtomFrozen != NULL) {202:         CudaSafeCall( cudaFree(m_d_isAtomFrozen) );
207:                 CudaSafeCall( cudaFree(m_d_isAtomFrozen) ); 
208:         } 
209: 203: 
210:         CudaSafeCall( cudaFree(m_d_nLargeRigidBody) );204:         CudaSafeCall( cudaFree(m_d_nLargeRigidBody) );
211:         CudaSafeCall( cudaFree(m_d_largeRigidIndices) );205:         CudaSafeCall( cudaFree(m_d_largeRigidIndices) );
212: 206: 
213:         CudaSafeCall( cudaFree(m_d_ljAddRep) );207:         CudaSafeCall( cudaFree(m_d_ljAddRep) );
214:         CudaSafeCall( cudaFree(m_d_ljAddAtt) );208:         CudaSafeCall( cudaFree(m_d_ljAddAtt) );
215: 209: 
216:         delete [] m_largeRigidIndices;210:         delete [] m_largeRigidIndices;
217: }211: }
218: 212: 


r32617/CUDAGMINOPTIM.cu 2017-05-26 12:30:25.904895302 +0100 r32616/CUDAGMINOPTIM.cu 2017-05-26 12:30:28.276926228 +0100
128: }128: }
129: 129: 
130: extern "C" void gminoptim_copy_forces(gpuContext gpu, double* gradient)130: extern "C" void gminoptim_copy_forces(gpuContext gpu, double* gradient)
131: {131: {
132:         dim3 blockDim;132:         dim3 blockDim;
133:         blockDim.x = 1024;133:         blockDim.x = 1024;
134:         dim3 gridDim;134:         dim3 gridDim;
135:         gridDim.x = ((gpu->sim.atoms)+ blockDim.x - 1)/blockDim.x;135:         gridDim.x = ((gpu->sim.atoms)+ blockDim.x - 1)/blockDim.x;
136:         get_frc<<<gridDim, blockDim>>>(gpu->pbForce->_pDevData, gradient);        136:         get_frc<<<gridDim, blockDim>>>(gpu->pbForce->_pDevData, gradient);        
137: }137: }
 138: 
 139: extern "C" void gminoptim_coords_cputogpu(gpuContext gpu, int *natoms, double *coords)
 140: {
 141:         double *d_coords;
 142:         cudaMalloc(&d_coords, 3 * (*natoms) * sizeof(double));
 143:         cudaMemcpy(d_coords, coords, 3 * (*natoms) * sizeof(double), cudaMemcpyHostToDevice);
 144:         gminoptim_copy_atoms(gpu, d_coords);
 145:         cudaFree(d_coords);
 146: }
 147: 
 148: extern "C" void gminoptim_ene_gputocpu(gpuContext gpu, double *energy)
 149: {
 150:         double *d_energy;
 151:         cudaMalloc(&d_energy, sizeof(double));
 152:         gminoptim_copy_ene(gpu, d_energy);
 153:         cudaMemcpy(energy, d_energy, sizeof(double), cudaMemcpyDeviceToHost);
 154:         cudaFree(d_energy);
 155: }
 156: 
 157: extern "C" void gminoptim_gradients_gputocpu(gpuContext gpu, int *natoms, double *gradients)
 158: {
 159:         double *d_gradients;
 160:         cudaMalloc(&d_gradients, 3 * (*natoms) * sizeof(double));
 161:         gminoptim_copy_forces(gpu, d_gradients);
 162:         cudaMemcpy(gradients, d_gradients, 3 * (*natoms) * sizeof(double), cudaMemcpyDeviceToHost);
 163:         cudaFree(d_gradients);
 164: }


r32617/gpu.cpp 2017-05-26 12:30:26.132898275 +0100 r32616/gpu.cpp 2017-05-26 12:30:28.508929254 +0100
827: }827: }
828: gpuCopyConstants();828: gpuCopyConstants();
829: }829: }
830: 830: 
831: // Function below added for CUDAGMIN/CUDAOPTIM831: // Function below added for CUDAGMIN/CUDAOPTIM
832: extern "C" void gminoptim_copy_const_()832: extern "C" void gminoptim_copy_const_()
833: {833: {
834:         gminoptim_copy_constants(gpu);834:         gminoptim_copy_constants(gpu);
835: }835: }
836: 836: 
 837: //Function below added for CUDAGMIN/CUDAOPTIM
 838: extern "C" void gminoptim_enegrad_cputogpu(int *natoms, double *coords, double *energy, double *gradients)
 839: {
 840:         gminoptim_coords_cputogpu(gpu, natoms, coords);
 841:         gminoptim_gpu_gb_ene();
 842:         gminoptim_ene_gputocpu(gpu, energy);
 843:         gminoptim_gradients_gputocpu(gpu, natoms, gradients);
 844: }
 845: 
837: // Function below added for CUDAGMIN/CUDAOPTIM846: // Function below added for CUDAGMIN/CUDAOPTIM
838: extern "C" void gminoptim_copy_crd(const double* atm_crd)847: extern "C" void gminoptim_copy_crd(const double* atm_crd)
839: {848: {
840:         gminoptim_copy_atoms(gpu, atm_crd);849:         gminoptim_copy_atoms(gpu, atm_crd);
841: }850: }
842: 851: 
843: extern "C" void gpu_upload_crd_(double atm_crd[][3])852: extern "C" void gpu_upload_crd_(double atm_crd[][3])
844: {853: {
845:         PRINTMETHOD("gpu_upload_crd");854:         PRINTMETHOD("gpu_upload_crd");
846:         if (gpu->bNeighborList && (gpu->pbImageIndex != NULL))855:         if (gpu->bNeighborList && (gpu->pbImageIndex != NULL))


r32617/keywords.f 2017-05-26 12:30:26.800906984 +0100 r32616/keywords.f 2017-05-26 12:30:29.176937963 +0100
1192:       MLQPROB=.FALSE.1192:       MLQPROB=.FALSE.
1193:       MLQDONE=.FALSE.1193:       MLQDONE=.FALSE.
1194:       MLQNORM=.FALSE.1194:       MLQNORM=.FALSE.
1195:       MLQLAMBDA=0.0D01195:       MLQLAMBDA=0.0D0
1196:       MLQSTART=11196:       MLQSTART=1
1197: 1197: 
1198:       LJADDT=.FALSE.1198:       LJADDT=.FALSE.
1199:       LJADD2T=.FALSE.1199:       LJADD2T=.FALSE.
1200:       LJADD3T=.FALSE.1200:       LJADD3T=.FALSE.
1201:       LJADD4T=.FALSE.1201:       LJADD4T=.FALSE.
1202:       NADDTARGET=01202:       NADDTARGET=1
1203:       REORDERADDT=.FALSE.1203:       REORDERADDT=.FALSE.
1204:       PYADDT=.FALSE.1204:       PYADDT=.FALSE.
1205:       PYADD2T=.FALSE.1205:       PYADD2T=.FALSE.
1206: 1206: 
1207:       DUMPMQT=.FALSE.1207:       DUMPMQT=.FALSE.
1208: ! jk669 SQNM keywords1208: ! jk669 SQNM keywords
1209:       SQNMT=.FALSE.1209:       SQNMT=.FALSE.
1210:       SQNM_HISTMAX=20 !defualt1210:       SQNM_HISTMAX=20 !defualt
1211:       SQNM_DEBUGT=.FALSE.1211:       SQNM_DEBUGT=.FALSE.
1212:       SQNM_DEBUGRUN=0 !defaul 0 means debug all runs1212:       SQNM_DEBUGRUN=0 !defaul 0 means debug all runs


r32617/lj_potential.cu 2017-05-26 12:30:26.572904013 +0100 r32616/lj_potential.cu 2017-05-26 12:30:28.948934992 +0100
265:                         r.x = myPositionX - refPositionX;265:                         r.x = myPositionX - refPositionX;
266:                         r.y = myPositionY - refPositionY;266:                         r.y = myPositionY - refPositionY;
267:                         r.z = myPositionZ - refPositionZ;267:                         r.z = myPositionZ - refPositionZ;
268: 268: 
269:                         // Distance squared. 269:                         // Distance squared. 
270:                         double rsq = r.x * r.x + r.y * r.y + r.z * r.z;270:                         double rsq = r.x * r.x + r.y * r.y + r.z * r.z;
271: 271: 
272:                         double attractive = m_d_ljAddAtt[mj2+mj1*nAddTarget];272:                         double attractive = m_d_ljAddAtt[mj2+mj1*nAddTarget];
273:                         double repulsive = m_d_ljAddRep[mj2+mj1*nAddTarget];273:                         double repulsive = m_d_ljAddRep[mj2+mj1*nAddTarget];
274: 274: 
275:                         // Ignore pair if distance is too small, this is usually the self-interaction of the atom. 275:                         d_energyContrib[tid] = 0.0;
276:                         double energyContrib = 0.0;276:                         d_gradientContrib[3*tid+0] = 0.0;
277:                         double gradContrib = 0.0;277:                         d_gradientContrib[3*tid+1] = 0.0;
 278:                         d_gradientContrib[3*tid+2] = 0.0;
278: 279: 
 280:                         // Ignore pair if distance is too small, this is usually the self-interaction of the atom. 
279:                         if (rsq >= 1.0e-6) {281:                         if (rsq >= 1.0e-6) {
280:                                 // Calculate 1/r**2, 1/r**6, 1/r**12. 282:                                 // Calculate 1/r**2, 1/r**6, 1/r**12. 
281:                                 double ir2  = 1.0 / rsq;283:                                 double ir2  = 1.0 / rsq;
282:                                 double ir6  = ir2*ir2*ir2;284:                                 double ir6  = ir2*ir2*ir2;
283:                                 double ir12 = ir6*ir6;285:                                 double ir12 = ir6*ir6;
284: 286: 
285:                                 // Calculate the energy. 287:                                 // Calculate the energy. 
286:                                 energyContrib = 2.0 * (ir12*repulsive-ir6*attractive);288:                                 d_energyContrib[tid] = 2.0 * (ir12*repulsive-ir6*attractive);
287:                                 // Calculate the gradient. 289:                                 // Calculate the gradient. 
288:                                 gradContrib = 4.0 * (12.0 * ir12 * repulsive -290:                                 double gradcontrib = 4.0 * (12.0 * ir12 * repulsive -
289:                                                 6.0 * ir6 * attractive) * ir2;291:                                                 6.0 * ir6 * attractive) * ir2;
 292:                                 d_gradientContrib[3*tid+0] = gradcontrib * r.x;
 293:                                 d_gradientContrib[3*tid+1] = gradcontrib * r.y;
 294:                                 d_gradientContrib[3*tid+2] = gradcontrib * r.z;
290:                         }295:                         }
291: 296: 
292:                         d_energyContrib[tid] = energyContrib; 
293:                         d_gradientContrib[3*tid+0] = gradContrib * r.x; 
294:                         d_gradientContrib[3*tid+1] = gradContrib * r.y; 
295:                         d_gradientContrib[3*tid+2] = gradContrib * r.z; 
296:  
297:                         tid += blockDim.x * gridDim.x;297:                         tid += blockDim.x * gridDim.x;
298:                 }298:                 }
299:         }299:         }
300: 300: 
301: }301: }


r32617/main.F 2017-05-26 12:30:27.168911862 +0100 r32616/main.F 2017-05-26 12:30:29.392940782 +0100
172:       IF (CALCQT) CALL SHINIT172:       IF (CALCQT) CALL SHINIT
173: 173: 
174:       ALLOCATE(FINISH(3*NATOMSALLOC))174:       ALLOCATE(FINISH(3*NATOMSALLOC))
175: 175: 
176:       INQUIRE(UNIT=1,OPENED=LOPEN)176:       INQUIRE(UNIT=1,OPENED=LOPEN)
177:       IF (LOPEN) THEN177:       IF (LOPEN) THEN
178:          WRITE(*,'(A,I2,A)') 'main> A ERROR *** Unit ', 1, ' is not free '178:          WRITE(*,'(A,I2,A)') 'main> A ERROR *** Unit ', 1, ' is not free '
179:          STOP179:          STOP
180:       ENDIF180:       ENDIF
181:       CALL KEYWORD181:       CALL KEYWORD
182:       IF ((NADDTARGET .EQ. 0) .AND. CUDAT) THEN 
183:          NADDTARGET = NATOMS 
184:       END IF 
185: !     write out the atom indices to check:182: !     write out the atom indices to check:
186: !        IF(TWISTT) THEN183: !        IF(TWISTT) THEN
187: !              WRITE(MYUNIT,*)'NWISTGROUPS dihedrals to be constrained..'184: !              WRITE(MYUNIT,*)'NWISTGROUPS dihedrals to be constrained..'
188: !              DO J1=1, NTWISTGROUPS185: !              DO J1=1, NTWISTGROUPS
189: !                WRITE(MYUNIT,*) 'Group J1  indices:'186: !                WRITE(MYUNIT,*) 'Group J1  indices:'
190: !                WRITE(MYUNIT,*) TWIST_ATOMS(1:4, J1)187: !                WRITE(MYUNIT,*) TWIST_ATOMS(1:4, J1)
191: !                WRITE(MYUNIT,*) 'Group J1  k value, ref.dihedral:'188: !                WRITE(MYUNIT,*) 'Group J1  k value, ref.dihedral:'
192: !                WRITE(MYUNIT,*) TWIST_K(J1), TWIST_THETA0(J1)189: !                WRITE(MYUNIT,*) TWIST_K(J1), TWIST_THETA0(J1)
193: !              END DO190: !              END DO
194: !        END IF191: !        END IF


r32617/modcudalbfgs.f90 2017-05-26 12:30:27.828920388 +0100 r32616/modcudalbfgs.f90 2017-05-26 12:30:30.060949490 +0100
 59:         REAL(KIND=C_DOUBLE), DIMENSION(C_NRIGIDBODY*3*3), INTENT(IN) :: C_IINVERSE ! RBF: inverse eigenvalues of the unweighted tensor of gyration 59:         REAL(KIND=C_DOUBLE), DIMENSION(C_NRIGIDBODY*3*3), INTENT(IN) :: C_IINVERSE ! RBF: inverse eigenvalues of the unweighted tensor of gyration
 60:  60: 
 61:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix 61:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix
 62:  62: 
 63:         REAL(KIND=C_DOUBLE), DIMENSION(N), INTENT(INOUT) :: C_XCOORDS ! Coordinates 63:         REAL(KIND=C_DOUBLE), DIMENSION(N), INTENT(INOUT) :: C_XCOORDS ! Coordinates
 64:  64: 
 65:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_ENERGY, & ! Energy 65:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_ENERGY, & ! Energy
 66:                                             C_RMS, & ! RMS force 66:                                             C_RMS, & ! RMS force
 67:                                             POTENTIALTIME ! Time taken in calculating potential 67:                                             POTENTIALTIME ! Time taken in calculating potential
 68:  68: 
 69:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_DEBUG, & ! If true, print debug info 69:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_DEBUG, & ! If true, print debug info.
 70:                                             C_CUDATIMET, & ! If true, print timing info 70:                                             C_CUDATIMET, & ! If true, print timing info.
 71:                                             C_ATOMRIGIDCOORDT, & ! If false, use rigid body coordinates 71:                                             C_ATOMRIGIDCOORDT, & ! If false, use rigid body coordinates
 72:                                             PROJECT, & ! If true, project out the components of the gradient along the eigenvector 72:                                             PROJECT, & ! If true, project out the components of the gradient along the eigenvector
 73:                                             C_FREEZE, & ! If true, freeze some specified atoms 73:                                             C_FREEZE, & ! If true, freeze some specified atoms
 74:                                             C_AACONVERGENCET ! If true, use more accurate method of calculating rigid RMS force 74:                                             C_AACONVERGENCET ! If true, use more accurate method of calculating rigid RMS force
 75:  75: 
 76:         LOGICAL(KIND=C_BOOL), DIMENSION(N/3), INTENT(IN) :: C_FROZEN ! Logical array specifying frozen atoms 76:         LOGICAL(KIND=C_BOOL), DIMENSION(N/3), INTENT(IN) :: C_FROZEN ! Logical array specifying frozen atoms
 77:  77: 
 78:         LOGICAL(KIND=C_BOOL), INTENT(OUT) :: C_MFLAG, & ! True if quench converged 78:         LOGICAL(KIND=C_BOOL), INTENT(OUT) :: C_MFLAG, & ! True if quench converged
 79:                                              C_COLDFUSION ! Set to true during minimization if cold fusion diagnosed 79:                                              C_COLDFUSION ! Set to true during minimization if cold fusion diagnosed
 80:  80: 
 81:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used 81:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used
 82:  82: 
 83:     END SUBROUTINE CUDA_LBFGS 83:     END SUBROUTINE CUDA_LBFGS
 84: END INTERFACE 84: END INTERFACE
 85:  85: 
 86: INTERFACE 86: INTERFACE
 87:     SUBROUTINE CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, &  87:     SUBROUTINE CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, & 
 88:                                     C_CUDATIMET, POTENTIALTIME) BIND(C,NAME="setup_potential_cputogpu") 88:                                      C_LJADDATT) BIND(C,NAME="gminoptim_enegrad_cputogpu")
 89:  89: 
 90:         IMPORT :: C_INT, C_DOUBLE, C_BOOL, C_CHAR 90:         IMPORT :: C_INT, C_DOUBLE
 91:  91: 
 92:         INTEGER(KIND=C_INT), INTENT(IN) :: NATOMS, & ! No. of atoms 92:         INTEGER(KIND=C_INT), INTENT(IN) :: NATOMS, & ! No. of atoms
 93:                                            C_NADDTARGET ! Target cluster size (addressability) 93:                                            C_NADDTARGET ! Target cluster size (addressability)
 94:  94: 
 95:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(IN) :: COORDS ! Atomic coordinates 95:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(IN) :: COORDS ! Atomic coordinates
 96:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix 96:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix
 97:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_TOTENERGY, & ! Total energy of the system 97:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_TOTENERGY ! Total energy of the system
 98:                                             POTENTIALTIME ! Time taken in calculating potential 
 99:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(OUT) :: C_GRADIENTS ! Gradient of the energy w.r.t. each atomic coordinate 98:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(OUT) :: C_GRADIENTS ! Gradient of the energy w.r.t. each atomic coordinate
100:  99: 
101:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_CUDATIMET ! If true, print timing info 
102:  
103:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used 
104:  
105:     END SUBROUTINE CUDA_ENEGRAD_CPUTOGPU100:     END SUBROUTINE CUDA_ENEGRAD_CPUTOGPU
106: END INTERFACE101: END INTERFACE
107: 102: 
108: CONTAINS103: CONTAINS
109: 104: 
110:     SUBROUTINE CUDA_LBFGS_WRAPPER(N, MUPDATE, XCOORDS, MFLAG, ENERGY, RMS, ITMAX, ITDONE, PROJECT)105:     SUBROUTINE CUDA_LBFGS_WRAPPER(N, MUPDATE, XCOORDS, MFLAG, ENERGY, RMS, ITMAX, ITDONE, PROJECT)
111: 106: 
112:         ! Variables passed as *arguments through this wrapper* (not common) with intent in for CUDA_LBFGS are converted directly107:         ! Variables passed as *arguments through this wrapper* (not common) with intent in for CUDA_LBFGS are converted directly
113:         INTEGER(KIND=C_INT) :: N, ITMAX, C_ITDONE, NCALLS, C_DEGFREEDOMS, C_NRIGIDBODY, C_MAXSITE, MUPDATE, C_NFREEZE, C_NADDTARGET108:         INTEGER(KIND=C_INT) :: N, ITMAX, C_ITDONE, NCALLS, C_DEGFREEDOMS, C_NRIGIDBODY, C_MAXSITE, MUPDATE, C_NFREEZE, C_NADDTARGET
114:         INTEGER(KIND=C_INT), DIMENSION(NRIGIDBODY) :: C_NSITEPERBODY109:         INTEGER(KIND=C_INT), DIMENSION(NRIGIDBODY) :: C_NSITEPERBODY
264:         KNOWG=.FALSE.259:         KNOWG=.FALSE.
265:         KNOWH=.FALSE.260:         KNOWH=.FALSE.
266: 261: 
267:     END SUBROUTINE CUDA_LBFGS_WRAPPER262:     END SUBROUTINE CUDA_LBFGS_WRAPPER
268: 263: 
269: 264: 
270:     SUBROUTINE CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, TOTENERGY, GRADIENTS)265:     SUBROUTINE CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, TOTENERGY, GRADIENTS)
271: 266: 
272:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET267:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET
273: 268: 
274:         REAL(KIND=C_DOUBLE) :: C_TOTENERGY, POTENTIALTIME269:         REAL(KIND=C_DOUBLE) :: C_TOTENERGY
275:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS) :: COORDS, C_GRADIENTS270:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS) :: COORDS, C_GRADIENTS
276:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT271:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT
277: 272: 
278:         LOGICAL(KIND=C_BOOL) :: C_CUDATIMET 
279:  
280:         CHARACTER(LEN=1, KIND=C_CHAR) :: C_CUDAPOT 
281:  
282:         INTEGER :: X, I, J273:         INTEGER :: X, I, J
283: 274: 
284:         DOUBLE PRECISION :: TOTENERGY275:         DOUBLE PRECISION :: TOTENERGY
285:         DOUBLE PRECISION, DIMENSION(3*NATOMS) :: GRADIENTS276:         DOUBLE PRECISION, DIMENSION(3*NATOMS) :: GRADIENTS
286: 277: 
287:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN278:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN
288:             DO J = 1,NADDTARGET279:             DO J = 1,NADDTARGET
289:                 DO I = 1,NADDTARGET280:                 DO I = 1,NADDTARGET
290:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)281:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)
291:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)282:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)
292:                 END DO283:                 END DO
293:             END DO284:             END DO
294:         ELSE285:         ELSE
295:             C_LJADDREP(:) = 1.0286:             C_LJADDREP(:) = 1.0
296:             C_LJADDATT(:) = 1.0287:             C_LJADDATT(:) = 1.0
297:         END IF288:         END IF
298: 289: 
299:         C_NADDTARGET = NADDTARGET290:         C_NADDTARGET = NADDTARGET
300:         C_CUDAPOT = CUDAPOT 
301:         C_CUDATIMET = CUDATIMET 
302: 291: 
303:         ! Calculates the energy and gradients on the GPU using the GB potential292:         ! Calculates the energy and gradients on the GPU using the GB potential
304:         CALL CUDA_ENEGRAD_CPUTOGPU(3*NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, & 293:         CALL CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT)
305:                                   C_CUDATIMET, POTENTIALTIME) 
306: 294: 
307:         TOTENERGY = DBLE(C_TOTENERGY)295:         TOTENERGY = DBLE(C_TOTENERGY)
308: 296: 
309:         DO X = 1,(3*NATOMS)297:         DO X = 1,(3*NATOMS)
310:             GRADIENTS(X) = DBLE(C_GRADIENTS(X))298:             GRADIENTS(X) = DBLE(C_GRADIENTS(X))
311:         END DO299:         END DO
312: 300: 
313:         ! FTIME and FCALL are updated in potential.f 
314:  
315:     END SUBROUTINE CUDA_ENEGRAD_WRAPPER301:     END SUBROUTINE CUDA_ENEGRAD_WRAPPER
316: 302: 
317: END MODULE MODCUDALBFGS303: END MODULE MODCUDALBFGS


r32617/modcudalbfgs.F90 2017-05-26 12:30:27.384914598 +0100 r32616/modcudalbfgs.F90 2017-05-26 12:30:29.616943700 +0100
 58:         REAL(KIND=C_DOUBLE), DIMENSION(C_NRIGIDBODY*3*3), INTENT(IN) :: C_IINVERSE ! RBF: inverse eigenvalues of the unweighted tensor of gyration 58:         REAL(KIND=C_DOUBLE), DIMENSION(C_NRIGIDBODY*3*3), INTENT(IN) :: C_IINVERSE ! RBF: inverse eigenvalues of the unweighted tensor of gyration
 59:  59: 
 60:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix 60:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix
 61:  61: 
 62:         REAL(KIND=C_DOUBLE), DIMENSION(N), INTENT(INOUT) :: C_XCOORDS ! Coordinates 62:         REAL(KIND=C_DOUBLE), DIMENSION(N), INTENT(INOUT) :: C_XCOORDS ! Coordinates
 63:  63: 
 64:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_ENERGY, & ! Energy 64:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_ENERGY, & ! Energy
 65:                                             C_RMS, & ! RMS force 65:                                             C_RMS, & ! RMS force
 66:                                             POTENTIALTIME ! Time taken in calculating potential - not used in GMIN 66:                                             POTENTIALTIME ! Time taken in calculating potential - not used in GMIN
 67:  67: 
 68:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_DEBUG, & ! If true, print debug info 68:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_DEBUG, & ! If true, print debug info.
 69:                                             C_CUDATIMET, & ! If true, print timing info 69:                                             C_CUDATIMET, & ! If true, print timing info. 
 70:                                             C_ATOMRIGIDCOORDT, & ! If false, use rigid body coordinates 70:                                             C_ATOMRIGIDCOORDT, & ! If false, use rigid body coordinates
 71:                                             C_FREEZE, & ! If true, freeze some specified atoms 71:                                             C_FREEZE, & ! If true, freeze some specified atoms
 72:                                             PROJECT, & ! PROJECT is OPTIM only, always false 72:                                             PROJECT, & ! PROJECT is OPTIM only, always false
 73:                                             C_AACONVERGENCET ! If true, use more accurate method of calculating rigid RMS force 73:                                             C_AACONVERGENCET ! If true, use more accurate method of calculating rigid RMS force
 74:  74: 
 75:         LOGICAL(KIND=C_BOOL), DIMENSION(N/3), INTENT(IN) :: C_FROZEN ! Logical array specifying frozen atoms 75:         LOGICAL(KIND=C_BOOL), DIMENSION(N/3), INTENT(IN) :: C_FROZEN ! Logical array specifying frozen atoms
 76:  76: 
 77:         LOGICAL(KIND=C_BOOL), INTENT(OUT) :: C_MFLAG, & ! True if quench converged 77:         LOGICAL(KIND=C_BOOL), INTENT(OUT) :: C_MFLAG, & ! True if quench converged
 78:                                              C_COLDFUSION ! Set to true during minimization if cold fusion diagnosed 78:                                              C_COLDFUSION ! Set to true during minimization if cold fusion diagnosed
 79:  79: 
 80:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used 80:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used
 81:  81: 
 82:     END SUBROUTINE CUDA_LBFGS 82:     END SUBROUTINE CUDA_LBFGS
 83: END INTERFACE 83: END INTERFACE
 84:  84: 
 85: INTERFACE 85: INTERFACE
 86:     SUBROUTINE CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, &  86:     SUBROUTINE CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, & 
 87:                                     C_CUDATIMET, POTENTIALTIME) BIND(C,NAME="setup_potential_cputogpu") 87:                                      C_LJADDATT) BIND(C,NAME="gminoptim_enegrad_cputogpu")
 88:  88: 
 89:         IMPORT :: C_INT, C_DOUBLE, C_BOOL, C_CHAR 89:         IMPORT :: C_INT, C_DOUBLE
 90:  90: 
 91:         INTEGER(KIND=C_INT), INTENT(IN) :: NATOMS, & ! No. of atoms 91:         INTEGER(KIND=C_INT), INTENT(IN) :: NATOMS, & ! No. of atoms
 92:                                            C_NADDTARGET ! Target cluster size (addressability) 92:                                            C_NADDTARGET ! Target cluster size (addressability)
 93:  93: 
 94:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(IN) :: COORDS ! Atomic coordinates 94:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(IN) :: COORDS ! Atomic coordinates
 95:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix 95:         REAL(KIND=C_DOUBLE), DIMENSION(C_NADDTARGET, C_NADDTARGET), INTENT(IN) ::C_LJADDREP, C_LJADDATT ! Repulsive/attractive epsilon matrix
 96:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_TOTENERGY, & ! Total energy of the system 96:         REAL(KIND=C_DOUBLE), INTENT(OUT) :: C_TOTENERGY ! Total energy of the system
 97:                                             POTENTIALTIME ! Time taken in calculating potential - not used in GMIN 
 98:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(OUT) :: C_GRADIENTS ! Gradient of the energy w.r.t. each atomic coordinate 97:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS), INTENT(OUT) :: C_GRADIENTS ! Gradient of the energy w.r.t. each atomic coordinate
 99:  98: 
100:         LOGICAL(KIND=C_BOOL), INTENT(IN) :: C_CUDATIMET ! If true, print timing info 
101:  
102:         CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: C_CUDAPOT ! Character specifying the CUDA potential to be used 
103:  
104:     END SUBROUTINE CUDA_ENEGRAD_CPUTOGPU     99:     END SUBROUTINE CUDA_ENEGRAD_CPUTOGPU    
105: END INTERFACE100: END INTERFACE
106: #endif /* DUMMY_CUDA */101: #endif /* DUMMY_CUDA */
107: 102: 
108: CONTAINS103: CONTAINS
109: 104: 
110:     SUBROUTINE CUDA_LBFGS_WRAPPER(N, MUPDATE, XCOORDS, EPS, MFLAG, ENERGY, ITMAX, ITDONE, RESET)105:     SUBROUTINE CUDA_LBFGS_WRAPPER(N, MUPDATE, XCOORDS, EPS, MFLAG, ENERGY, ITMAX, ITDONE, RESET)
111: 106: 
112:         ! Variables passed as *arguments through this wrapper* (not common) with intent in for CUDA_LBFGS are converted directly107:         ! Variables passed as *arguments through this wrapper* (not common) with intent in for CUDA_LBFGS are converted directly
113: #ifndef DUMMY_CUDA108: #ifndef DUMMY_CUDA
277:         END IF272:         END IF
278: 273: 
279: #endif /* DUMMY_CUDA */274: #endif /* DUMMY_CUDA */
280: 275: 
281:     END SUBROUTINE CUDA_LBFGS_WRAPPER276:     END SUBROUTINE CUDA_LBFGS_WRAPPER
282: 277: 
283:     SUBROUTINE CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, TOTENERGY, GRADIENTS)278:     SUBROUTINE CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, TOTENERGY, GRADIENTS)
284: 279: 
285:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET280:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET
286: 281: 
287:         REAL(KIND=C_DOUBLE) :: C_TOTENERGY, POTENTIALTIME282:         REAL(KIND=C_DOUBLE) :: C_TOTENERGY
288: #ifndef DUMMY_CUDA 
289:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS) :: COORDS, C_GRADIENTS283:         REAL(KIND=C_DOUBLE), DIMENSION(3*NATOMS) :: COORDS, C_GRADIENTS
 284: #ifndef DUMMY_CUDA
290:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT285:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT
291: #else ifdef DUMMY_CUDA286: #else /* ifdef DUMMY_CUDA */
292:         REAL(KIND=C_DOUBLE), DIMENSION(1) :: COORDS, C_GRADIENTS 
293:         REAL(KIND=C_DOUBLE), DIMENSION(1, 1) :: C_LJADDREP, C_LJADDATT287:         REAL(KIND=C_DOUBLE), DIMENSION(1, 1) :: C_LJADDREP, C_LJADDATT
294: #endif288: #endif
295:  
296:         LOGICAL(KIND=C_BOOL) :: C_CUDATIMET 
297:  
298:         CHARACTER(LEN=1, KIND=C_CHAR) :: C_CUDAPOT 
299:  
300:         INTEGER :: X, I, J289:         INTEGER :: X, I, J
301: 290: 
302:         DOUBLE PRECISION :: TOTENERGY291:         DOUBLE PRECISION :: TOTENERGY
303:         DOUBLE PRECISION, DIMENSION(3*NATOMS) :: GRADIENTS292:         DOUBLE PRECISION, DIMENSION(3*NATOMS) :: GRADIENTS
304: 293: 
305: #ifndef DUMMY_CUDA294: #ifndef DUMMY_CUDA
306: 295: 
307:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN296:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN
308:             DO J = 1,NADDTARGET297:             DO J = 1,NADDTARGET
309:                 DO I = 1,NADDTARGET298:                 DO I = 1,NADDTARGET
310:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)299:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)
311:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)300:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)
312:                 END DO301:                 END DO
313:             END DO302:             END DO
314:         ELSE303:         ELSE
315:             C_LJADDREP(:) = 1.0304:             C_LJADDREP(:) = 1.0
316:             C_LJADDATT(:) = 1.0305:             C_LJADDATT(:) = 1.0
317:         END IF306:         END IF
318: 307: 
319:         C_NADDTARGET = NADDTARGET308:         C_NADDTARGET = NADDTARGET
320:         C_CUDAPOT = CUDAPOT 
321:         C_CUDATIMET = CUDATIMET 
322: 309: 
323:         ! Calculates the energy and gradients on the GPU using the GB potential310:         ! Calculates the energy and gradients on the GPU using the GB potential
324:         CALL CUDA_ENEGRAD_CPUTOGPU(3*NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, & 311:         CALL CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_TOTENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT)
325:                                   C_CUDATIMET, POTENTIALTIME) 
326: 312: 
327:         TOTENERGY = DBLE(C_TOTENERGY)313:         TOTENERGY = DBLE(C_TOTENERGY)
328: 314:         
329:         DO X = 1,(3*NATOMS)315:         DO X = 1,(3*NATOMS)
330:             GRADIENTS(X) = DBLE(C_GRADIENTS(X))316:             GRADIENTS(X) = DBLE(C_GRADIENTS(X))
331:         END DO317:         END DO
332: #endif318: #endif
333:     END SUBROUTINE CUDA_ENEGRAD_WRAPPER319:     END SUBROUTINE CUDA_ENEGRAD_WRAPPER
334: 320: 
335:     SUBROUTINE CUDA_NUMERICAL_HESS(NATOMS, COORDS, HESSIAN, DELTA)321:     SUBROUTINE CUDA_NUMERICAL_HESS(NATOMS, COORDS, HESSIAN, DELTA)
336:         IMPLICIT NONE322:         IMPLICIT NONE
337: 323: 
338:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET324:         INTEGER(KIND=C_INT) :: NATOMS, C_NADDTARGET
339:         REAL(KIND=C_DOUBLE) :: C_ENERGY, POTENTIALTIME325:         REAL(KIND=C_DOUBLE) :: C_ENERGY
340:         REAL(KIND=C_DOUBLE) :: COORDS(3*NATOMS), C_GRADIENTS(3*NATOMS)326:         REAL(KIND=C_DOUBLE) :: COORDS(3*NATOMS), C_GRADIENTS(3*NATOMS)
341: #ifndef DUMMY_CUDA327: #ifndef DUMMY_CUDA
342:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT328:         REAL(KIND=C_DOUBLE), DIMENSION(NADDTARGET*NADDTARGET) :: C_LJADDREP, C_LJADDATT
343: #else ifdef DUMMY_CUDA329: #else /* ifdef DUMMY_CUDA */
344:         REAL(KIND=C_DOUBLE), DIMENSION(1, 1) :: C_LJADDREP, C_LJADDATT330:         REAL(KIND=C_DOUBLE), DIMENSION(1, 1) :: C_LJADDREP, C_LJADDATT
345: #endif331: #endif
346:  
347:         LOGICAL(KIND=C_BOOL) :: C_CUDATIMET 
348:  
349:         CHARACTER(LEN=1, KIND=C_CHAR) :: C_CUDAPOT 
350:  
351:         DOUBLE PRECISION    :: HESSIAN(3*NATOMS, 3*NATOMS)332:         DOUBLE PRECISION    :: HESSIAN(3*NATOMS, 3*NATOMS)
352:         DOUBLE PRECISION    :: DELTA333:         DOUBLE PRECISION    :: DELTA
353:         DOUBLE PRECISION    :: GRAD_PLUS(3*NATOMS), GRAD_MINUS(3*NATOMS)334:         DOUBLE PRECISION    :: GRAD_PLUS(3*NATOMS), GRAD_MINUS(3*NATOMS)
354:         INTEGER             :: I, J335:         INTEGER             :: I, J
355: 336: 
356: #ifndef DUMMY_CUDA337: #ifndef DUMMY_CUDA
357:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN338:         IF (ALLOCATED(LJADDREP) .AND. ALLOCATED(LJADDATT) .AND. LJADD3T) THEN
358:             DO J = 1,NADDTARGET339:             DO J = 1,NADDTARGET
359:                 DO I = 1,NADDTARGET340:                 DO I = 1,NADDTARGET
360:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)341:                     C_LJADDREP((J - 1)*NADDTARGET + I) = LJADDREP(I,J)
361:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)342:                     C_LJADDATT((J - 1)*NADDTARGET + I) = LJADDATT(I,J)
362:                 END DO343:                 END DO
363:             END DO344:             END DO
364:         ELSE345:         ELSE
365:             C_LJADDREP(:) = 1.0346:             C_LJADDREP(:) = 1.0
366:             C_LJADDATT(:) = 1.0347:             C_LJADDATT(:) = 1.0
367:         END IF348:         END IF
368: 349: 
369:         C_NADDTARGET = NADDTARGET350:         C_NADDTARGET = NADDTARGET
370:         C_CUDAPOT = CUDAPOT 
371:         C_CUDATIMET = CUDATIMET 
372: 351: 
373:         DO I = 1, 3*NATOMS352:         DO I = 1, 3*NATOMS
374:             ! Plus353:             ! Plus
375:             COORDS(I) = COORDS(I) + DELTA354:             COORDS(I) = COORDS(I) + DELTA
376:             CALL CUDA_ENEGRAD_CPUTOGPU(3*NATOMS, COORDS, C_ENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, & 355:             CALL CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_ENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT)
377:                                       C_CUDATIMET, POTENTIALTIME) 
378:             GRAD_PLUS(:) = DBLE(C_GRADIENTS(:))356:             GRAD_PLUS(:) = DBLE(C_GRADIENTS(:))
379:             ! Minus357:             ! Minus
380:             COORDS(I) = COORDS(I) - 2.0D0 * DELTA358:             COORDS(I) = COORDS(I) - 2.0D0 * DELTA
381:             CALL CUDA_ENEGRAD_CPUTOGPU(3*NATOMS, COORDS, C_ENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT, C_CUDAPOT, & 359:             CALL CUDA_ENEGRAD_CPUTOGPU(NATOMS, COORDS, C_ENERGY, C_GRADIENTS, C_NADDTARGET, C_LJADDREP, C_LJADDATT)
382:                                       C_CUDATIMET, POTENTIALTIME) 
383:             GRAD_MINUS(:) = DBLE(C_GRADIENTS(:))360:             GRAD_MINUS(:) = DBLE(C_GRADIENTS(:))
384:             ! Reset coords361:             ! Reset coords
385:             COORDS(I) = COORDS(I) + DELTA362:             COORDS(I) = COORDS(I) + DELTA
386:             ! Calculate hessian363:             ! Calculate hessian
387:             HESSIAN(I, :) = (GRAD_PLUS(:) - GRAD_MINUS(:)) / (2.0D0 * DELTA)364:             HESSIAN(I, :) = (GRAD_PLUS(:) - GRAD_MINUS(:)) / (2.0D0 * DELTA)
388:         END DO365:         END DO
389: #endif366: #endif
390:     END SUBROUTINE CUDA_NUMERICAL_HESS367:     END SUBROUTINE CUDA_NUMERICAL_HESS
391: 368: 
392: END MODULE MODCUDALBFGS369: END MODULE MODCUDALBFGS


r32617/potential.f 2017-05-26 12:30:28.056923361 +0100 r32616/potential.f 2017-05-26 12:30:30.292952515 +0100
1082: !!                      ENDIF1082: !!                      ENDIF
1083: !!                   ENDDO1083: !!                   ENDDO
1084: !!                ENDDO1084: !!                ENDDO
1085: !                 STOP1085: !                 STOP
1086: 1086: 
1087:             IF (PTEST) THEN1087:             IF (PTEST) THEN
1088:                WRITE(*,20) ' potential> Energy for last cycle=',ENERGY1088:                WRITE(*,20) ' potential> Energy for last cycle=',ENERGY
1089:                WRITE(ESTRING,20) ' potential> Energy for last cycle=',ENERGY1089:                WRITE(ESTRING,20) ' potential> Energy for last cycle=',ENERGY
1090:             ENDIF1090:             ENDIF
1091:          ELSE IF (LJADD3T) THEN1091:          ELSE IF (LJADD3T) THEN
1092:             IF (CUDAT) THEN1092:             CALL LJADD3(NATOMS, COORDS, VNEW, ENERGY, GTEST, SSTEST)
1093:                ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU 
1094:                CALL CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, ENERGY, GRADATOMS) 
1095:                VNEW(1:3*NATOMS) = GRADATOMS(:) 
1096:             ELSE 
1097:                CALL LJADD3(NATOMS, COORDS, VNEW, ENERGY, GTEST, SSTEST) 
1098:             END IF 
1099:             IF (PTEST) THEN1093:             IF (PTEST) THEN
1100:                WRITE(*,20) ' potential> Energy for last cycle=',ENERGY1094:                WRITE(*,20) ' potential> Energy for last cycle=',ENERGY
1101:                WRITE(ESTRING,20) ' potential> Energy for last cycle=',ENERGY1095:                WRITE(ESTRING,20) ' potential> Energy for last cycle=',ENERGY
1102:             ENDIF1096:             ENDIF
1103:          ELSE IF (LJADD2T) THEN1097:          ELSE IF (LJADD2T) THEN
1104:             CALL LJADD2(NATOMS, COORDS, VNEW, ENERGY, GTEST, SSTEST)1098:             CALL LJADD2(NATOMS, COORDS, VNEW, ENERGY, GTEST, SSTEST)
1105: !                DIFF=1.0D-41099: !                DIFF=1.0D-4
1106: !                PRINT*,'analytic and numerical gradients:'1100: !                PRINT*,'analytic and numerical gradients:'
1107: !                IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(3*NATOMS,3*NATOMS))1101: !                IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(3*NATOMS,3*NATOMS))
1108: !                CALL LJADD2(NATOMS, COORDS, VNEW, ENERGY, .TRUE., .TRUE.)1102: !                CALL LJADD2(NATOMS, COORDS, VNEW, ENERGY, .TRUE., .TRUE.)


r32617/potential.f90 2017-05-26 12:30:27.608917521 +0100 r32616/potential.f90 2017-05-26 12:30:29.840946622 +0100
182:          RETURN182:          RETURN
183:       END IF183:       END IF
184:    ELSE IF (LJADDT) THEN184:    ELSE IF (LJADDT) THEN
185: 185: 
186:       IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN186:       IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN
187:           XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS)187:           XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS)
188:           CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS)188:           CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS)
189:       ENDIF189:       ENDIF
190: 190: 
191:       IF (LJADD3T) THEN191:       IF (LJADD3T) THEN
192:          IF (CUDAT) THEN192:          CALL LJADD3(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT)
193: ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU 
194:             CALL CUDA_ENEGRAD_WRAPPER(NATOMS, X, EREAL, GRADATOMS) 
195:          ELSE 
196:             CALL LJADD3(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT) 
197:          END IF 
198:       ELSEIF (LJADD2T) THEN193:       ELSEIF (LJADD2T) THEN
199:          CALL LJADD2(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT)194:          CALL LJADD2(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT)
200:       ELSE195:       ELSE
201:          CALL LJADD(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT)196:          CALL LJADD(NATOMS, X, GRADATOMS, EREAL, GRADT, SECT)
202:       ENDIF197:       ENDIF
203:       GRAD(1:3*NATOMS)=GRADATOMS(:)198:       GRAD(1:3*NATOMS)=GRADATOMS(:)
204: 199: 
205:       IF ( RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN200:       IF ( RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN
206: 201: 
207:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0202:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0
492: !         GRAD(J1)=(EPLUS-EMINUS)/(2.0D0*DIFF)487: !         GRAD(J1)=(EPLUS-EMINUS)/(2.0D0*DIFF)
493: !      END DO488: !      END DO
494: !      WRITE(6,*)'SUMM DIFF ', SUMMDIFF489: !      WRITE(6,*)'SUMM DIFF ', SUMMDIFF
495: 490: 
496: ! khs26> AMBER12 energy and gradient call491: ! khs26> AMBER12 energy and gradient call
497:    ELSE IF (AMBER12T) THEN492:    ELSE IF (AMBER12T) THEN
498: ! If coordinates are atomistic (i.e. no rigid bodies) just call the energy and gradient493: ! If coordinates are atomistic (i.e. no rigid bodies) just call the energy and gradient
499:       IF (ATOMRIGIDCOORDT) THEN494:       IF (ATOMRIGIDCOORDT) THEN
500:          IF (CUDAT) THEN495:          IF (CUDAT) THEN
501: ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU496: ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU
 497: ! Calls to CUDA_ENEGRAD_WRAPPER are for debugging/testing/printing purposes
502:             CALL CUDA_ENEGRAD_WRAPPER(NATOMS, X, EREAL, GRADATOMS)498:             CALL CUDA_ENEGRAD_WRAPPER(NATOMS, X, EREAL, GRADATOMS)
503:             GRAD(1:3*NATOMS)=GRADATOMS(:)499:             GRAD(1:3*NATOMS)=GRADATOMS(:)
504:          ELSE500:          ELSE
505:             CALL AMBER12_ENERGY_AND_GRADIENT(NATOMS, X, EREAL, GRADATOMS, AMBER12_ENERGY_DECOMP)501:             CALL AMBER12_ENERGY_AND_GRADIENT(NATOMS, X, EREAL, GRADATOMS, AMBER12_ENERGY_DECOMP)
506:             GRAD(1:3*NATOMS)=GRADATOMS(:)502:             GRAD(1:3*NATOMS)=GRADATOMS(:)
507:          END IF503:          END IF
508: ! If the coordinates include rigid bodies, transform them back to being atomistic first504: ! If the coordinates include rigid bodies, transform them back to being atomistic first
509:       ELSE505:       ELSE
510:          XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS)506:          XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS)
511:          CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, XCOORDS, XRIGIDCOORDS)507:          CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, XCOORDS, XRIGIDCOORDS)
522:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)518:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)
523:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0519:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0
524:          GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)520:          GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)
525:          GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0521:          GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0
526:       END IF522:       END IF
527: !copied from OPTIM, rbody part not tested523: !copied from OPTIM, rbody part not tested
528:       IF (SECT) THEN524:       IF (SECT) THEN
529:          IF (ALLOCATED(HESS)) DEALLOCATE(HESS)525:          IF (ALLOCATED(HESS)) DEALLOCATE(HESS)
530:          ALLOCATE(HESS(3*NATOMS, 3*NATOMS))526:          ALLOCATE(HESS(3*NATOMS, 3*NATOMS))
531:          IF (CUDAT) THEN527:          IF (CUDAT) THEN
532:             CALL CUDA_NUMERICAL_HESS(NATOMS, X, HESS, DELTA=1.0D-4)528:             CALL CUDA_NUMERICAL_HESS(NATOMS, X, HESS, DELTA=1.0D-6)
533:          ELSE529:          ELSE
534:             CALL AMBER12_NUM_HESS(NATOMS,X, DELTA=1.0D-4, HESSIAN=HESS(:, :))530:             CALL AMBER12_NUM_HESS(NATOMS,X, DELTA=1.0D-5, HESSIAN=HESS(:, :))
535:          END IF531:          END IF
536:          IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT) ) THEN532:          IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT) ) THEN
537:             CALL TRANSFORMHESSIAN(HESS, GRADATOMS, XRIGIDCOORDS,XRIGIDHESS, RBAANORMALMODET)533:             CALL TRANSFORMHESSIAN(HESS, GRADATOMS, XRIGIDCOORDS,XRIGIDHESS, RBAANORMALMODET)
538:             HESS(DEGFREEDOMS+1:3*NATOMS,:) = 0.0D0534:             HESS(DEGFREEDOMS+1:3*NATOMS,:) = 0.0D0
539:             HESS(:,DEGFREEDOMS+1:3*NATOMS) = 0.0D0535:             HESS(:,DEGFREEDOMS+1:3*NATOMS) = 0.0D0
540:             HESS(1:DEGFREEDOMS,1:DEGFREEDOMS) = XRIGIDHESS(1:DEGFREEDOMS,1:DEGFREEDOMS)536:             HESS(1:DEGFREEDOMS,1:DEGFREEDOMS) = XRIGIDHESS(1:DEGFREEDOMS,1:DEGFREEDOMS)
541:          END IF537:          END IF
542:       END IF538:       END IF
543: ! AMBER 9 Energy and gradient calls539: ! AMBER 9 Energy and gradient calls
544:    ELSE IF (AMBERT) THEN540:    ELSE IF (AMBERT) THEN
1352:    ELSE1348:    ELSE
1353: ! Otherwise, assume Lennard-Jones type particles.1349: ! Otherwise, assume Lennard-Jones type particles.
1354: 1350: 
1355: ! RAD must be called before the routine that calculates the potential or LBFGS1351: ! RAD must be called before the routine that calculates the potential or LBFGS
1356: ! will get confused even if EVAP is set .TRUE. correctly.1352: ! will get confused even if EVAP is set .TRUE. correctly.
1357:       CALL RAD(X, GRAD, EREAL, GRADT)1353:       CALL RAD(X, GRAD, EREAL, GRADT)
1358:       IF (EVAPREJECT) RETURN1354:       IF (EVAPREJECT) RETURN
1359:       IF (CUTT) THEN1355:       IF (CUTT) THEN
1360:          CALL LJCUT(X, GRAD, EREAL, GRADT, SECT)1356:          CALL LJCUT(X, GRAD, EREAL, GRADT, SECT)
1361:       ELSE1357:       ELSE
1362:          IF (CUDAT) THEN1358:          CALL LJ(X, GRAD, EREAL, GRADT, SECT)
1363: ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU 
1364:             CALL CUDA_ENEGRAD_WRAPPER(NATOMS, X, EREAL, GRADATOMS) 
1365:             GRAD(1:3*NATOMS)=GRADATOMS(:) 
1366:          ELSE 
1367:             CALL LJ(X, GRAD, EREAL, GRADT, SECT) 
1368:          END IF 
1369:          IF (AXTELL) CALL AXT(NATOMS, X, GRAD, EREAL, GRADT, ZSTAR)1359:          IF (AXTELL) CALL AXT(NATOMS, X, GRAD, EREAL, GRADT, ZSTAR)
1370:       END IF1360:       END IF
1371:    END IF1361:    END IF
1372: !1362: !
1373: !  --------------- End of possible potentials - now add fields if required ------------------------------1363: !  --------------- End of possible potentials - now add fields if required ------------------------------
1374: !1364: !
1375: 1365: 
1376:    IF (RIGIDCONTOURT) THEN1366:    IF (RIGIDCONTOURT) THEN
1377: !       sf344> sanity check, since this keyword works only for a pair of rigid particles1367: !       sf344> sanity check, since this keyword works only for a pair of rigid particles
1378:       IF (NATOMS.NE.4) THEN1368:       IF (NATOMS.NE.4) THEN


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0