|
CP2K 2.4 (Revision 12889)
|
00001 !-----------------------------------------------------------------------------! 00002 ! CP2K: A general program to perform molecular dynamics simulations ! 00003 ! Copyright (C) 2000 - 2013 CP2K developers group ! 00004 !-----------------------------------------------------------------------------! 00005 00006 ! ***************************************************************************** 00014 MODULE cp_subsys_types 00015 00016 USE atomic_kind_list_types, ONLY: atomic_kind_list_release,& 00017 atomic_kind_list_retain,& 00018 atomic_kind_list_type 00019 USE cell_types, ONLY: cell_type,& 00020 real_to_scaled,& 00021 scaled_to_real 00022 USE colvar_types, ONLY: colvar_p_type,& 00023 colvar_release 00024 USE cp_para_env, ONLY: cp_para_env_release,& 00025 cp_para_env_retain 00026 USE cp_para_types, ONLY: cp_para_env_type 00027 USE distribution_1d_types, ONLY: distribution_1d_release,& 00028 distribution_1d_retain,& 00029 distribution_1d_type 00030 USE f77_blas 00031 USE kinds, ONLY: dp 00032 USE mol_kind_new_list_types, ONLY: mol_kind_new_list_release,& 00033 mol_kind_new_list_retain,& 00034 mol_kind_new_list_type 00035 USE mol_new_list_types, ONLY: mol_new_list_release,& 00036 mol_new_list_retain,& 00037 mol_new_list_type 00038 USE molecule_types_new, ONLY: global_constraint_type 00039 USE multipole_types, ONLY: multipole_type,& 00040 release_multipole_type,& 00041 retain_multipole_type 00042 USE particle_list_types, ONLY: particle_list_release,& 00043 particle_list_retain,& 00044 particle_list_type 00045 #include "cp_common_uses.h" 00046 00047 IMPLICIT NONE 00048 PRIVATE 00049 00050 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_subsys_types' 00051 INTEGER, PRIVATE, SAVE :: last_fragment_id = 0 00052 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .FALSE. 00053 00054 PUBLIC :: cp_subsys_type,& 00055 cp_subsys_p_type 00056 00057 PUBLIC :: cp_subsys_create,& 00058 cp_subsys_retain,& 00059 cp_subsys_release,& 00060 cp_subsys_get,& 00061 cp_subsys_set,& 00062 pack_subsys_particles,& 00063 unpack_subsys_particles 00064 00065 ! ***************************************************************************** 00080 TYPE cp_subsys_type 00081 INTEGER :: ref_count, id_nr 00082 TYPE (atomic_kind_list_type), POINTER :: atomic_kinds 00083 TYPE (particle_list_type), POINTER :: particles 00084 TYPE (particle_list_type), POINTER :: shell_particles 00085 TYPE (particle_list_type), POINTER :: core_particles 00086 TYPE (distribution_1d_type), POINTER :: local_particles 00087 TYPE (cp_para_env_type), POINTER :: para_env 00088 ! New molecules kinds 00089 TYPE (mol_new_list_type), POINTER :: molecules_new 00090 TYPE (mol_kind_new_list_type), POINTER :: molecule_kinds_new 00091 TYPE (distribution_1d_type), POINTER :: local_molecules_new 00092 ! Definitions of the collective variables 00093 TYPE (colvar_p_type), DIMENSION(:), POINTER :: colvar_p 00094 ! Intermolecular constraints 00095 TYPE (global_constraint_type), POINTER :: gci 00096 ! Multipoles 00097 TYPE (multipole_type), POINTER :: multipoles 00098 END TYPE cp_subsys_type 00099 00100 ! ***************************************************************************** 00108 TYPE cp_subsys_p_type 00109 TYPE(cp_subsys_type), POINTER :: subsys 00110 END TYPE cp_subsys_p_type 00111 00112 CONTAINS 00113 00114 ! ***************************************************************************** 00124 SUBROUTINE cp_subsys_create(subsys, para_env, error) 00125 TYPE(cp_subsys_type), POINTER :: subsys 00126 TYPE(cp_para_env_type), POINTER :: para_env 00127 TYPE(cp_error_type), INTENT(inout) :: error 00128 00129 CHARACTER(len=*), PARAMETER :: routineN = 'cp_subsys_create', 00130 routineP = moduleN//':'//routineN 00131 00132 INTEGER :: stat 00133 LOGICAL :: failure 00134 00135 failure=.FALSE. 00136 00137 ALLOCATE(subsys, stat=stat) 00138 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00139 IF (.NOT. failure) THEN 00140 last_fragment_id=last_fragment_id+1 00141 subsys%id_nr=last_fragment_id 00142 subsys%ref_count=1 00143 CALL cp_para_env_retain(para_env,error=error) 00144 subsys%para_env => para_env 00145 NULLIFY(subsys%atomic_kinds, subsys%particles, & 00146 subsys%shell_particles , subsys%core_particles, & 00147 subsys%local_particles, subsys%molecules_new,& 00148 subsys%molecule_kinds_new, subsys%local_molecules_new,& 00149 subsys%gci, subsys%multipoles) 00150 NULLIFY(subsys%colvar_p) 00151 END IF 00152 END SUBROUTINE cp_subsys_create 00153 00154 ! ***************************************************************************** 00163 SUBROUTINE cp_subsys_retain(subsys, error) 00164 TYPE(cp_subsys_type), POINTER :: subsys 00165 TYPE(cp_error_type), INTENT(inout) :: error 00166 00167 CHARACTER(len=*), PARAMETER :: routineN = 'cp_subsys_retain', 00168 routineP = moduleN//':'//routineN 00169 00170 LOGICAL :: failure 00171 00172 failure=.FALSE. 00173 00174 CPPrecondition(ASSOCIATED(subsys),cp_failure_level,routineP,error,failure) 00175 IF (.NOT. failure) THEN 00176 CPPreconditionNoFail(subsys%ref_count>0,cp_failure_level,routineP,error) 00177 subsys%ref_count=subsys%ref_count+1 00178 END IF 00179 END SUBROUTINE cp_subsys_retain 00180 00181 ! ***************************************************************************** 00190 SUBROUTINE cp_subsys_release(subsys, error) 00191 TYPE(cp_subsys_type), POINTER :: subsys 00192 TYPE(cp_error_type), INTENT(inout) :: error 00193 00194 CHARACTER(len=*), PARAMETER :: routineN = 'cp_subsys_release', 00195 routineP = moduleN//':'//routineN 00196 00197 INTEGER :: i, j, stat 00198 LOGICAL :: failure 00199 00200 failure=.FALSE. 00201 00202 IF (ASSOCIATED(subsys)) THEN 00203 CPPreconditionNoFail(subsys%ref_count>0,cp_failure_level,routineP,error) 00204 subsys%ref_count=subsys%ref_count-1 00205 IF (subsys%ref_count==0) THEN 00206 CALL atomic_kind_list_release(subsys%atomic_kinds,error=error) 00207 CALL particle_list_release(subsys%particles, error=error) 00208 CALL particle_list_release(subsys%shell_particles, error=error) 00209 CALL particle_list_release(subsys%core_particles, error=error) 00210 CALL distribution_1d_release(subsys%local_particles, error=error) 00211 CALL mol_kind_new_list_release(subsys%molecule_kinds_new, error=error) 00212 CALL mol_new_list_release(subsys%molecules_new, error=error) 00213 CALL distribution_1d_release(subsys%local_molecules_new,error=error) 00214 CALL cp_para_env_release(subsys%para_env, error=error) 00215 ! Multipoles 00216 IF(ASSOCIATED(subsys%multipoles)) THEN 00217 CALL release_multipole_type(subsys%multipoles, error) 00218 END IF 00219 ! Colvar info 00220 IF(ASSOCIATED(subsys%colvar_p)) THEN 00221 DO i=1,SIZE(subsys%colvar_p) 00222 IF (ASSOCIATED(subsys%colvar_p(i)%colvar)) THEN 00223 CALL colvar_release(subsys%colvar_p(i)%colvar,error) 00224 END IF 00225 ENDDO 00226 DEALLOCATE(subsys%colvar_p,STAT=stat) 00227 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00228 END IF 00229 ! Intermolecular constraints 00230 IF (ASSOCIATED(subsys%gci)) THEN 00231 ! List of constraints 00232 IF (ASSOCIATED(subsys%gci%colv_list)) THEN 00233 DO j = 1, SIZE(subsys%gci%colv_list) 00234 DEALLOCATE (subsys%gci%colv_list(j)%i_atoms,STAT=stat) 00235 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00236 END DO 00237 DEALLOCATE (subsys%gci%colv_list,STAT=stat) 00238 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00239 END IF 00240 IF (ASSOCIATED(subsys%gci%g3x3_list)) THEN 00241 DEALLOCATE (subsys%gci%g3x3_list,STAT=stat) 00242 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00243 END IF 00244 IF (ASSOCIATED(subsys%gci%g4x6_list)) THEN 00245 DEALLOCATE (subsys%gci%g4x6_list,STAT=stat) 00246 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00247 END IF 00248 ! Local information 00249 IF (ASSOCIATED(subsys%gci%lcolv)) THEN 00250 DO j = 1, SIZE(subsys%gci%lcolv) 00251 CALL colvar_release(subsys%gci%lcolv(j)%colvar,error=error) 00252 CALL colvar_release(subsys%gci%lcolv(j)%colvar_old,error=error) 00253 NULLIFY(subsys%gci%lcolv(j)%colvar) 00254 NULLIFY(subsys%gci%lcolv(j)%colvar_old) 00255 END DO 00256 DEALLOCATE (subsys%gci%lcolv,STAT=stat) 00257 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00258 ENDIF 00259 IF (ASSOCIATED(subsys%gci%lg3x3)) THEN 00260 DEALLOCATE (subsys%gci%lg3x3,STAT=stat) 00261 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00262 ENDIF 00263 IF (ASSOCIATED(subsys%gci%lg4x6)) THEN 00264 DEALLOCATE (subsys%gci%lg4x6,STAT=stat) 00265 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00266 ENDIF 00267 IF (ASSOCIATED(subsys%gci%fixd_list)) THEN 00268 DEALLOCATE (subsys%gci%fixd_list,STAT=stat) 00269 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00270 ENDIF 00271 DEALLOCATE (subsys%gci,STAT=stat) 00272 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00273 ENDIF 00274 DEALLOCATE(subsys, stat=stat) 00275 CPPostconditionNoFail(stat==0,cp_warning_level,routineP,error) 00276 END IF 00277 END IF 00278 END SUBROUTINE cp_subsys_release 00279 00280 ! ***************************************************************************** 00291 SUBROUTINE cp_subsys_set(subsys, atomic_kinds, particles, local_particles,& 00292 molecules_new, molecule_kinds_new, local_molecules_new, para_env,& 00293 colvar_p, shell_particles, core_particles, gci, multipoles, error) 00294 TYPE(cp_subsys_type), POINTER :: subsys 00295 TYPE(atomic_kind_list_type), OPTIONAL, 00296 POINTER :: atomic_kinds 00297 TYPE(particle_list_type), OPTIONAL, 00298 POINTER :: particles 00299 TYPE(distribution_1d_type), OPTIONAL, 00300 POINTER :: local_particles 00301 TYPE(mol_new_list_type), OPTIONAL, 00302 POINTER :: molecules_new 00303 TYPE(mol_kind_new_list_type), OPTIONAL, 00304 POINTER :: molecule_kinds_new 00305 TYPE(distribution_1d_type), OPTIONAL, 00306 POINTER :: local_molecules_new 00307 TYPE(cp_para_env_type), OPTIONAL, 00308 POINTER :: para_env 00309 TYPE(colvar_p_type), DIMENSION(:), 00310 OPTIONAL, POINTER :: colvar_p 00311 TYPE(particle_list_type), OPTIONAL, 00312 POINTER :: shell_particles, 00313 core_particles 00314 TYPE(global_constraint_type), OPTIONAL, 00315 POINTER :: gci 00316 TYPE(multipole_type), OPTIONAL, POINTER :: multipoles 00317 TYPE(cp_error_type), INTENT(inout) :: error 00318 00319 CHARACTER(len=*), PARAMETER :: routineN = 'cp_subsys_set', 00320 routineP = moduleN//':'//routineN 00321 00322 LOGICAL :: failure 00323 00324 failure=.FALSE. 00325 00326 CPPrecondition(ASSOCIATED(subsys),cp_failure_level,routineP,error,failure) 00327 CPPrecondition(subsys%ref_count>0,cp_failure_level,routineP,error,failure) 00328 IF (.NOT. failure) THEN 00329 IF (PRESENT(multipoles)) THEN 00330 CALL retain_multipole_type(multipoles, error) 00331 CALL release_multipole_type(subsys%multipoles, error) 00332 subsys%multipoles => multipoles 00333 END IF 00334 IF (PRESENT(atomic_kinds)) THEN 00335 CALL atomic_kind_list_retain(atomic_kinds,error=error) 00336 CALL atomic_kind_list_release(subsys%atomic_kinds, error=error) 00337 subsys%atomic_kinds => atomic_kinds 00338 END IF 00339 IF (PRESENT(particles)) THEN 00340 CALL particle_list_retain(particles, error=error) 00341 CALL particle_list_release(subsys%particles, error=error) 00342 subsys%particles => particles 00343 END IF 00344 IF (PRESENT(local_particles)) THEN 00345 CALL distribution_1d_retain(local_particles,error=error) 00346 CALL distribution_1d_release(subsys%local_particles,error=error) 00347 subsys%local_particles => local_particles 00348 END IF 00349 IF (PRESENT(local_molecules_new)) THEN 00350 CALL distribution_1d_retain(local_molecules_new,error=error) 00351 CALL distribution_1d_release(subsys%local_molecules_new,error=error) 00352 subsys%local_molecules_new => local_molecules_new 00353 END IF 00354 IF (PRESENT(molecule_kinds_new)) THEN 00355 CALL mol_kind_new_list_retain(molecule_kinds_new, error=error) 00356 CALL mol_kind_new_list_release(subsys%molecule_kinds_new, error=error) 00357 subsys%molecule_kinds_new => molecule_kinds_new 00358 END IF 00359 IF (PRESENT(molecules_new)) THEN 00360 CALL mol_new_list_retain(molecules_new, error=error) 00361 CALL mol_new_list_release(subsys%molecules_new, error=error) 00362 subsys%molecules_new => molecules_new 00363 END IF 00364 IF (PRESENT(para_env)) THEN 00365 CALL cp_para_env_retain(para_env, error=error) 00366 CALL cp_para_env_release(subsys%para_env, error=error) 00367 subsys%para_env => para_env 00368 END IF 00369 IF (PRESENT(colvar_p)) THEN 00370 CPPrecondition(.NOT.ASSOCIATED(subsys%colvar_p),cp_failure_level,routineP,error,failure) 00371 subsys%colvar_p=>colvar_p 00372 ENDIF 00373 IF (PRESENT(shell_particles)) THEN 00374 IF(ASSOCIATED(shell_particles)) THEN 00375 CALL particle_list_retain(shell_particles, error=error) 00376 CALL particle_list_release(subsys%shell_particles, error=error) 00377 subsys%shell_particles => shell_particles 00378 END IF 00379 END IF 00380 IF (PRESENT(core_particles)) THEN 00381 IF(ASSOCIATED(core_particles)) THEN 00382 CALL particle_list_retain(core_particles, error=error) 00383 CALL particle_list_release(subsys%core_particles, error=error) 00384 subsys%core_particles => core_particles 00385 END IF 00386 END IF 00387 IF (PRESENT(gci)) THEN 00388 CPPrecondition(.NOT.ASSOCIATED(subsys%gci),cp_failure_level,routineP,error,failure) 00389 subsys%gci => gci 00390 ENDIF 00391 END IF 00392 END SUBROUTINE cp_subsys_set 00393 00394 ! ***************************************************************************** 00406 SUBROUTINE cp_subsys_get(subsys, id_nr, ref_count, atomic_kinds, particles,& 00407 local_particles, molecules_new, molecule_kinds_new,& 00408 local_molecules_new, para_env, colvar_p,& 00409 shell_particles, core_particles, gci, multipoles,& 00410 natom, nparticle, ncore, nshell, error) 00411 TYPE(cp_subsys_type), POINTER :: subsys 00412 INTEGER, INTENT(out), OPTIONAL :: id_nr, ref_count 00413 TYPE(atomic_kind_list_type), OPTIONAL, 00414 POINTER :: atomic_kinds 00415 TYPE(particle_list_type), OPTIONAL, 00416 POINTER :: particles 00417 TYPE(distribution_1d_type), OPTIONAL, 00418 POINTER :: local_particles 00419 TYPE(mol_new_list_type), OPTIONAL, 00420 POINTER :: molecules_new 00421 TYPE(mol_kind_new_list_type), OPTIONAL, 00422 POINTER :: molecule_kinds_new 00423 TYPE(distribution_1d_type), OPTIONAL, 00424 POINTER :: local_molecules_new 00425 TYPE(cp_para_env_type), OPTIONAL, 00426 POINTER :: para_env 00427 TYPE(colvar_p_type), DIMENSION(:), 00428 OPTIONAL, POINTER :: colvar_p 00429 TYPE(particle_list_type), OPTIONAL, 00430 POINTER :: shell_particles, 00431 core_particles 00432 TYPE(global_constraint_type), OPTIONAL, 00433 POINTER :: gci 00434 TYPE(multipole_type), OPTIONAL, POINTER :: multipoles 00435 INTEGER, INTENT(out), OPTIONAL :: natom, nparticle, ncore, 00436 nshell 00437 TYPE(cp_error_type), INTENT(inout) :: error 00438 00439 CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_subsys_get', 00440 routineP = moduleN//':'//routineN 00441 00442 INTEGER :: n_atom, n_core, n_shell 00443 LOGICAL :: failure 00444 00445 failure = .FALSE. 00446 n_atom = 0 00447 n_core = 0 00448 n_shell = 0 00449 00450 CPPrecondition(ASSOCIATED(subsys),cp_failure_level,routineP,error,failure) 00451 CPPrecondition(subsys%ref_count>0,cp_failure_level,routineP,error,failure) 00452 00453 IF (.NOT.failure) THEN 00454 IF (PRESENT(id_nr)) id_nr = subsys%id_nr 00455 IF (PRESENT(ref_count)) ref_count = subsys%ref_count 00456 IF (PRESENT(atomic_kinds)) atomic_kinds => subsys%atomic_kinds 00457 IF (PRESENT(particles)) particles => subsys%particles 00458 IF (PRESENT(local_particles)) local_particles => subsys%local_particles 00459 IF (PRESENT(molecules_new)) molecules_new => subsys%molecules_new 00460 IF (PRESENT(molecule_kinds_new)) molecule_kinds_new => subsys%molecule_kinds_new 00461 IF (PRESENT(local_molecules_new)) local_molecules_new => subsys%local_molecules_new 00462 IF (PRESENT(para_env)) para_env => subsys%para_env 00463 IF (PRESENT(colvar_p)) colvar_p => subsys%colvar_p 00464 IF (PRESENT(shell_particles)) shell_particles => subsys%shell_particles 00465 IF (PRESENT(core_particles)) core_particles => subsys%core_particles 00466 IF (PRESENT(gci)) gci => subsys%gci 00467 IF (PRESENT(multipoles)) multipoles => subsys%multipoles 00468 IF (PRESENT(natom).OR.PRESENT(nparticle).OR.PRESENT(nshell)) THEN 00469 ! An atomic particle set should be present in each subsystem at the moment 00470 CPPrecondition(ASSOCIATED(subsys%particles),cp_failure_level,routineP,error,failure) 00471 n_atom = subsys%particles%n_els 00472 ! Check if we have other kinds of particles in this subsystem 00473 IF (ASSOCIATED(subsys%shell_particles)) THEN 00474 n_shell = subsys%shell_particles%n_els 00475 CPPrecondition(ASSOCIATED(subsys%core_particles),cp_failure_level,routineP,error,failure) 00476 n_core = subsys%core_particles%n_els 00477 ! The same number of shell and core particles is assumed 00478 CPPrecondition((n_core == n_shell),cp_failure_level,routineP,error,failure) 00479 ELSE IF (ASSOCIATED(subsys%core_particles)) THEN 00480 ! This case should not occur at the moment 00481 CPPrecondition(ASSOCIATED(subsys%shell_particles),cp_failure_level,routineP,error,failure) 00482 ELSE 00483 n_core = 0 00484 n_shell = 0 00485 END IF 00486 IF (PRESENT(natom)) natom = n_atom 00487 IF (PRESENT(nparticle)) nparticle = n_atom + n_shell 00488 IF (PRESENT(ncore)) ncore = n_core 00489 IF (PRESENT(nshell)) nshell = n_shell 00490 END IF 00491 END IF 00492 00493 END SUBROUTINE cp_subsys_get 00494 00495 ! ***************************************************************************** 00502 SUBROUTINE pack_subsys_particles(subsys,f,r,s,v,fscale,cell,error) 00503 00504 TYPE(cp_subsys_type), POINTER :: subsys 00505 REAL(KIND=dp), DIMENSION(:), 00506 INTENT(OUT), OPTIONAL :: f, r, s, v 00507 REAL(KIND=dp), INTENT(IN), OPTIONAL :: fscale 00508 TYPE(cell_type), OPTIONAL, POINTER :: cell 00509 TYPE(cp_error_type), INTENT(INOUT) :: error 00510 00511 CHARACTER(LEN=*), PARAMETER :: routineN = 'pack_subsys_particles', 00512 routineP = moduleN//':'//routineN 00513 00514 INTEGER :: i, iatom, j, k, natom, 00515 nparticle, nsize, shell_index 00516 LOGICAL :: failure 00517 REAL(KIND=dp), DIMENSION(3) :: rs 00518 TYPE(particle_list_type), POINTER :: core_particles, particles, 00519 shell_particles 00520 00521 failure = .FALSE. 00522 CPPrecondition(ASSOCIATED(subsys),cp_failure_level,routineP,error,failure) 00523 00524 IF (PRESENT(s)) THEN 00525 CPPrecondition(PRESENT(cell),cp_failure_level,routineP,error,failure) 00526 END IF 00527 00528 NULLIFY (core_particles) 00529 NULLIFY (particles) 00530 NULLIFY (shell_particles) 00531 00532 CALL cp_subsys_get(subsys,& 00533 core_particles=core_particles,& 00534 natom=natom,& 00535 nparticle=nparticle,& 00536 particles=particles,& 00537 shell_particles=shell_particles,& 00538 error=error) 00539 00540 nsize = 3*nparticle 00541 00542 ! Pack forces 00543 00544 IF (PRESENT(f)) THEN 00545 CPPrecondition((SIZE(f) >= nsize),cp_failure_level,routineP,error,failure) 00546 j = 0 00547 DO iatom=1,natom 00548 shell_index = particles%els(iatom)%shell_index 00549 IF (shell_index == 0) THEN 00550 DO i=1,3 00551 j = j + 1 00552 f(j) = particles%els(iatom)%f(i) 00553 END DO 00554 ELSE 00555 DO i=1,3 00556 j = j + 1 00557 f(j) = core_particles%els(shell_index)%f(i) 00558 END DO 00559 k = 3*(natom + shell_index - 1) 00560 DO i=1,3 00561 f(k+i) = shell_particles%els(shell_index)%f(i) 00562 END DO 00563 END IF 00564 END DO 00565 IF (PRESENT(fscale)) f(1:nsize) = fscale*f(1:nsize) 00566 END IF 00567 00568 ! Pack coordinates 00569 00570 IF (PRESENT(r)) THEN 00571 CPPrecondition((SIZE(r) >= nsize),cp_failure_level,routineP,error,failure) 00572 j = 0 00573 DO iatom=1,natom 00574 shell_index = particles%els(iatom)%shell_index 00575 IF (shell_index == 0) THEN 00576 DO i=1,3 00577 j = j + 1 00578 r(j) = particles%els(iatom)%r(i) 00579 END DO 00580 ELSE 00581 DO i=1,3 00582 j = j + 1 00583 r(j) = core_particles%els(shell_index)%r(i) 00584 END DO 00585 k = 3*(natom + shell_index - 1) 00586 DO i=1,3 00587 r(k+i) = shell_particles%els(shell_index)%r(i) 00588 END DO 00589 END IF 00590 END DO 00591 END IF 00592 00593 ! Pack as scaled coordinates 00594 00595 IF (PRESENT(s)) THEN 00596 CPPrecondition((SIZE(s) >= nsize),cp_failure_level,routineP,error,failure) 00597 CPPrecondition(PRESENT(cell),cp_failure_level,routineP,error,failure) 00598 j = 0 00599 DO iatom=1,natom 00600 shell_index = particles%els(iatom)%shell_index 00601 IF (shell_index == 0) THEN 00602 CALL real_to_scaled(rs,particles%els(iatom)%r,cell) 00603 DO i=1,3 00604 j = j + 1 00605 s(j) = rs(i) 00606 END DO 00607 ELSE 00608 CALL real_to_scaled(rs,core_particles%els(shell_index)%r,cell) 00609 DO i=1,3 00610 j = j + 1 00611 s(j) = rs(i) 00612 END DO 00613 CALL real_to_scaled(rs,shell_particles%els(shell_index)%r,cell) 00614 k = 3*(natom + shell_index - 1) 00615 DO i=1,3 00616 s(k+i) = rs(i) 00617 END DO 00618 END IF 00619 END DO 00620 END IF 00621 00622 ! Pack velocities 00623 00624 IF (PRESENT(v)) THEN 00625 CPPrecondition((SIZE(v) >= nsize),cp_failure_level,routineP,error,failure) 00626 j = 0 00627 DO iatom=1,natom 00628 shell_index = particles%els(iatom)%shell_index 00629 IF (shell_index == 0) THEN 00630 DO i=1,3 00631 j = j + 1 00632 v(j) = particles%els(iatom)%v(i) 00633 END DO 00634 ELSE 00635 DO i=1,3 00636 j = j + 1 00637 v(j) = core_particles%els(shell_index)%v(i) 00638 END DO 00639 k = 3*(natom + shell_index - 1) 00640 DO i=1,3 00641 v(k+i) = shell_particles%els(shell_index)%v(i) 00642 END DO 00643 END IF 00644 END DO 00645 END IF 00646 00647 END SUBROUTINE pack_subsys_particles 00648 00649 ! ***************************************************************************** 00655 SUBROUTINE unpack_subsys_particles(subsys,f,r,s,v,fscale,cell,error) 00656 00657 TYPE(cp_subsys_type), POINTER :: subsys 00658 REAL(KIND=dp), DIMENSION(:), 00659 INTENT(IN), OPTIONAL :: f, r, s, v 00660 REAL(KIND=dp), INTENT(IN), OPTIONAL :: fscale 00661 TYPE(cell_type), OPTIONAL, POINTER :: cell 00662 TYPE(cp_error_type), INTENT(INOUT) :: error 00663 00664 CHARACTER(LEN=*), PARAMETER :: routineN = 'unpack_subsys_particles', 00665 routineP = moduleN//':'//routineN 00666 00667 INTEGER :: i, iatom, j, k, natom, 00668 nparticle, nsize, shell_index 00669 LOGICAL :: failure 00670 REAL(KIND=dp) :: fc, fs, mass, my_fscale 00671 REAL(KIND=dp), DIMENSION(3) :: rs 00672 TYPE(particle_list_type), POINTER :: core_particles, particles, 00673 shell_particles 00674 00675 failure = .FALSE. 00676 CPPrecondition(ASSOCIATED(subsys),cp_failure_level,routineP,error,failure) 00677 00678 NULLIFY (core_particles) 00679 NULLIFY (particles) 00680 NULLIFY (shell_particles) 00681 00682 CALL cp_subsys_get(subsys,& 00683 core_particles=core_particles,& 00684 natom=natom,& 00685 nparticle=nparticle,& 00686 particles=particles,& 00687 shell_particles=shell_particles,& 00688 error=error) 00689 00690 nsize = 3*nparticle 00691 00692 ! Unpack forces 00693 00694 IF (PRESENT(f)) THEN 00695 CPPrecondition((SIZE(f) >= nsize),cp_failure_level,routineP,error,failure) 00696 IF (PRESENT(fscale)) THEN 00697 my_fscale = fscale 00698 ELSE 00699 my_fscale = 1.0_dp 00700 END IF 00701 j = 0 00702 DO iatom=1,natom 00703 shell_index = particles%els(iatom)%shell_index 00704 IF (shell_index == 0) THEN 00705 DO i=1,3 00706 j = j + 1 00707 particles%els(iatom)%f(i) = my_fscale*f(j) 00708 END DO 00709 ELSE 00710 DO i=1,3 00711 j = j + 1 00712 core_particles%els(shell_index)%f(i) = my_fscale*f(j) 00713 END DO 00714 k = 3*(natom + shell_index - 1) 00715 DO i=1,3 00716 shell_particles%els(shell_index)%f(i) = my_fscale*f(k+i) 00717 END DO 00718 END IF 00719 END DO 00720 END IF 00721 00722 ! Unpack coordinates 00723 00724 IF (PRESENT(r)) THEN 00725 CPPrecondition((SIZE(r) >= nsize),cp_failure_level,routineP,error,failure) 00726 j = 0 00727 DO iatom=1,natom 00728 shell_index = particles%els(iatom)%shell_index 00729 IF (shell_index == 0) THEN 00730 DO i=1,3 00731 j = j + 1 00732 particles%els(iatom)%r(i) = r(j) 00733 END DO 00734 ELSE 00735 DO i=1,3 00736 j = j + 1 00737 core_particles%els(shell_index)%r(i) = r(j) 00738 END DO 00739 k = 3*(natom + shell_index - 1) 00740 DO i=1,3 00741 shell_particles%els(shell_index)%r(i) = r(k+i) 00742 END DO 00743 ! Update atomic position due to core and shell motion 00744 mass = particles%els(iatom)%atomic_kind%mass 00745 fc = core_particles%els(shell_index)%atomic_kind%shell%mass_core/mass 00746 fs = shell_particles%els(shell_index)%atomic_kind%shell%mass_shell/mass 00747 particles%els(iatom)%r(1:3) = fc*core_particles%els(shell_index)%r(1:3) +& 00748 fs*shell_particles%els(shell_index)%r(1:3) 00749 END IF 00750 END DO 00751 END IF 00752 00753 ! Unpack scaled coordinates 00754 00755 IF (PRESENT(s)) THEN 00756 CPPrecondition((SIZE(s) >= nsize),cp_failure_level,routineP,error,failure) 00757 CPPrecondition(PRESENT(cell),cp_failure_level,routineP,error,failure) 00758 j = 0 00759 DO iatom=1,natom 00760 shell_index = particles%els(iatom)%shell_index 00761 IF (shell_index == 0) THEN 00762 DO i=1,3 00763 j = j + 1 00764 rs(i) = s(j) 00765 END DO 00766 CALL scaled_to_real(particles%els(iatom)%r,rs,cell) 00767 ELSE 00768 DO i=1,3 00769 j = j + 1 00770 rs(i) = s(j) 00771 END DO 00772 CALL scaled_to_real(core_particles%els(shell_index)%r,rs,cell) 00773 k = 3*(natom + shell_index - 1) 00774 DO i=1,3 00775 rs(i) = s(k+i) 00776 END DO 00777 CALL scaled_to_real(shell_particles%els(shell_index)%r,rs,cell) 00778 ! Update atomic position due to core and shell motion 00779 mass = particles%els(iatom)%atomic_kind%mass 00780 fc = core_particles%els(shell_index)%atomic_kind%shell%mass_core/mass 00781 fs = shell_particles%els(shell_index)%atomic_kind%shell%mass_shell/mass 00782 particles%els(iatom)%r(1:3) = fc*core_particles%els(shell_index)%r(1:3) +& 00783 fs*shell_particles%els(shell_index)%r(1:3) 00784 END IF 00785 END DO 00786 END IF 00787 00788 ! Unpack velocities 00789 00790 IF (PRESENT(v)) THEN 00791 CPPrecondition((SIZE(v) >= nsize),cp_failure_level,routineP,error,failure) 00792 j = 0 00793 DO iatom=1,natom 00794 shell_index = particles%els(iatom)%shell_index 00795 IF (shell_index == 0) THEN 00796 DO i=1,3 00797 j = j + 1 00798 particles%els(iatom)%v(i) = v(j) 00799 END DO 00800 ELSE 00801 DO i=1,3 00802 j = j + 1 00803 core_particles%els(shell_index)%v(i) = v(j) 00804 END DO 00805 k = 3*(natom + shell_index - 1) 00806 DO i=1,3 00807 shell_particles%els(shell_index)%v(i) = v(k+i) 00808 END DO 00809 ! Update atomic velocity due to core and shell motion 00810 mass = particles%els(iatom)%atomic_kind%mass 00811 fc = core_particles%els(shell_index)%atomic_kind%shell%mass_core/mass 00812 fs = shell_particles%els(shell_index)%atomic_kind%shell%mass_shell/mass 00813 particles%els(iatom)%v(1:3) = fc*core_particles%els(shell_index)%v(1:3) +& 00814 fs*shell_particles%els(shell_index)%v(1:3) 00815 END IF 00816 END DO 00817 END IF 00818 00819 END SUBROUTINE unpack_subsys_particles 00820 00821 END MODULE cp_subsys_types
1.7.3