TABLE OF CONTENTS
ABINIT/pred_delocint [ Functions ]
NAME
pred_delocint
FUNCTION
Ionmov predictors (10) BFGS with delocalized internal coordinates IONMOV 10: Given a starting point xred that is a vector of length 3*(natom-1) (reduced nuclei coordinates), and unit cell parameters (acell and rprimd) the Broyden-Fletcher-Goldfarb-Shanno minimization is performed on the total energy function, using its gradient (atomic forces and stresses) while the optimization of unit cell parameters is only performed if optcell/=0. The convergence requirement on the atomic forces, 'tolmxf', allows an early exit. Otherwise no more than 'ntime' steps are performed. Returned quantities are xred, and eventually acell and rprimd (new ones!). Could see Numerical Recipes (Fortran), 1986, page 307. Implements the delocalized internal coordinate scheme of Andzelm et al. in CPL .335. 321 (2001) \ and Baker et al. JCP .105. 192 (1996) B matrix is derivative of delocalized internals wrt cartesian coordinates U matrix is eigenvectors of G = B*B^{T} S matrix is eigenvectors of F = B^{T}B
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, JCC, SE) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
ab_mover <type(abimover)> : Datatype with all the information needed by the preditor itime : Index of the present iteration ntime : Maximal number of iterations ionmov : (10 or 11) Specific kind of BFGS zDEBUG : if true print some debugging information
OUTPUT
SIDE EFFECTS
hist <type(abihist)> : History of positions,forces acell, rprimd, stresses
PARENTS
mover
CHILDREN
brdene,calc_prim_int,deloc2xcart,fcart2fred,fred2fdeloc,hessupdt hist2var,make_prim_internals,metric,var2hist,wrtout,xcart2deloc xcart2xred,xfh_recover_deloc,xfpack_f2vout,xfpack_vin2x,xfpack_x2vin xred2xcart
SOURCE
65 #if defined HAVE_CONFIG_H 66 #include "config.h" 67 #endif 68 69 #include "abi_common.h" 70 71 subroutine pred_delocint(ab_mover,ab_xfh,forstr,hist,ionmov,itime,zDEBUG,iexit) 72 73 use defs_basis 74 use m_profiling_abi 75 use m_abimover 76 use m_abihist 77 78 use m_geometry, only : fcart2fred, xcart2xred, xred2xcart, metric 79 use m_bfgs, only : hessinit, hessupdt, brdene 80 81 !This section has been created automatically by the script Abilint (TD). 82 !Do not modify the following lines by hand. 83 #undef ABI_FUNC 84 #define ABI_FUNC 'pred_delocint' 85 use interfaces_14_hidewrite 86 use interfaces_45_geomoptim, except_this_one => pred_delocint 87 !End of the abilint section 88 89 implicit none 90 91 !Arguments ------------------------------------ 92 !scalars 93 type(abimover),intent(inout) :: ab_mover 94 type(ab_xfh_type),intent(inout) :: ab_xfh 95 type(abihist),intent(inout) :: hist 96 type(abiforstr),intent(in) :: forstr 97 integer,intent(in) :: itime 98 integer,intent(in) :: ionmov 99 integer,intent(in) :: iexit 100 logical,intent(in) :: zDEBUG 101 102 !Local variables------------------------------- 103 !scalars 104 integer :: ndim,cycl_main 105 integer :: ihist_prev,ii,jj,kk 106 real(dp),save :: ucvol0 107 real(dp) :: ucvol 108 real(dp) :: etotal,etotal_prev 109 logical :: DEBUG=.TRUE. 110 integer,save :: icenter,irshift ! DELOCINT indexes 111 integer,save :: nshell,ndeloc ! DELOCINT number of 112 character(len=500) :: message 113 114 !arrays 115 real(dp),allocatable,save :: hessin(:,:),vin(:),vin_prev(:) 116 real(dp),allocatable,save :: vout(:),vout_prev(:) 117 real(dp),save :: acell0(3) ! Initial acell 118 real(dp),save :: rprimd0(3,3) ! Initial rprimd 119 real(dp),allocatable :: prim_int(:) 120 real(dp),allocatable,save :: u_matrix(:,:) ! DELOCINT this may need to be added to type inside ab_mover 121 real(dp) :: acell(3) 122 real(dp) :: rprimd(3,3) 123 real(dp) :: gprimd(3,3) 124 real(dp) :: gmet(3,3) 125 real(dp) :: rmet(3,3) 126 real(dp) :: residual(3,ab_mover%natom) 127 !real(dp) :: residual_corrected(3,ab_mover%natom) 128 real(dp) :: xred(3,ab_mover%natom),xcart(3,ab_mover%natom) 129 real(dp) :: strten(6) 130 real(dp) :: deloc_force(3*(ab_mover%natom-1)) 131 real(dp) :: deloc_int(3*(ab_mover%natom-1)) 132 real(dp) :: bt_inv_matrix(3*(ab_mover%natom-1),3*ab_mover%natom) 133 134 !*************************************************************************** 135 !Beginning of executable session 136 !*************************************************************************** 137 138 if(iexit/=0)then 139 if (allocated(vin)) then 140 ABI_DEALLOCATE(vin) 141 end if 142 if (allocated(vout)) then 143 ABI_DEALLOCATE(vout) 144 end if 145 if (allocated(vin_prev)) then 146 ABI_DEALLOCATE(vin_prev) 147 end if 148 if (allocated(vout_prev)) then 149 ABI_DEALLOCATE(vout_prev) 150 end if 151 if (allocated(hessin)) then 152 ABI_DEALLOCATE(hessin) 153 end if 154 if (allocated(u_matrix)) then 155 ABI_DEALLOCATE(u_matrix) 156 end if 157 return 158 end if 159 160 !write(std_out,*) 'delocint 01' 161 !########################################################## 162 !### 01. Debugging and Verbose 163 164 if(DEBUG)then 165 write(std_out,'(a,3a,38a,39a)') ch10,('-',kk=1,3),& 166 & 'Debugging and Verbose for pred_deloint',('-',kk=1,39) 167 write(std_out,*) 'ionmov: ',ionmov 168 write(std_out,*) 'itime: ',itime 169 end if 170 171 !write(std_out,*) 'delocint 02' 172 !########################################################## 173 !### 02. Compute the dimension of vectors (ndim) 174 175 !With internal we have 1 coordinate less 176 ndeloc = 3*(ab_mover%natom-1) 177 ndim=ndeloc 178 deloc_int(:)=zero 179 deloc_force(:)=zero 180 if(ab_mover%optcell==1 .or.& 181 & ab_mover%optcell==4 .or.& 182 & ab_mover%optcell==5 .or.& 183 & ab_mover%optcell==6) ndim=ndim+1 184 if(ab_mover%optcell==2 .or.& 185 & ab_mover%optcell==3) ndim=ndim+6 186 if(ab_mover%optcell==7 .or.& 187 & ab_mover%optcell==8 .or.& 188 & ab_mover%optcell==9) ndim=ndim+3 189 190 if(DEBUG) write(std_out,*) 'Dimension of vin, vout and hessian (ndim): ',ndim 191 192 !write(std_out,*) 'delocint 03' 193 !########################################################## 194 !### 03. Allocate the vectors vin, vout and hessian matrix 195 196 !Notice thqt vin, vout, etc could be allocated 197 !From a previous dataset with a different ndim 198 if(itime==1)then 199 if (allocated(vin)) then 200 ABI_DEALLOCATE(vin) 201 end if 202 if (allocated(vout)) then 203 ABI_DEALLOCATE(vout) 204 end if 205 if (allocated(vin_prev)) then 206 ABI_DEALLOCATE(vin_prev) 207 end if 208 if (allocated(vout_prev)) then 209 ABI_DEALLOCATE(vout_prev) 210 end if 211 if (allocated(hessin)) then 212 ABI_DEALLOCATE(hessin) 213 end if 214 ABI_ALLOCATE(vin,(ndim)) 215 ABI_ALLOCATE(vout,(ndim)) 216 ABI_ALLOCATE(vin_prev,(ndim)) 217 ABI_ALLOCATE(vout_prev,(ndim)) 218 ABI_ALLOCATE(hessin,(ndim,ndim)) 219 ! DELOCINT 220 ! Allocate all the variables of deloc, 221 ! Needed for compilers such as Pathscale that fails if 222 ! unallocated variables are passed as arguments 223 ! TODO: 224 ! all of the following should go into an init routine in m_abimover 225 226 nshell=3 227 ab_mover%deloc%nrshift=(2*nshell+1)**3 228 icenter = nshell*(2*nshell+1)**2 + nshell*(2*nshell+1) + nshell + 1 229 230 ABI_ALLOCATE(ab_mover%deloc%rshift,(3,ab_mover%deloc%nrshift)) 231 irshift=0 232 do ii=-nshell,nshell 233 do jj=-nshell,nshell 234 do kk=-nshell,nshell 235 irshift=irshift+1 236 ab_mover%deloc%rshift(:,irshift) = (/dble(ii),dble(jj),dble(kk)/) 237 end do 238 end do 239 end do 240 241 end if 242 243 244 !write(std_out,*) 'delocint 04' 245 !########################################################## 246 !### 04. Obtain the present values from the history 247 248 call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 249 call xred2xcart(ab_mover%natom,rprimd,xcart,xred) 250 251 strten(:)=hist%strten(:,hist%ihist) 252 etotal =hist%etot(hist%ihist) 253 254 !Fill the residual with forces (No preconditioning) 255 !Or the preconditioned forces 256 if (ab_mover%goprecon==0)then 257 call fcart2fred(hist%fcart(:,:,hist%ihist),residual,rprimd,ab_mover%natom) 258 else 259 residual(:,:)= forstr%fred(:,:) 260 end if 261 262 if(zDEBUG)then 263 write (std_out,*) 'residual:' 264 do kk=1,ab_mover%natom 265 write (std_out,*) residual(:,kk) 266 end do 267 write (std_out,*) 'strten:' 268 write (std_out,*) strten(1:3),ch10,strten(4:6) 269 write (std_out,*) 'etotal:' 270 write (std_out,*) etotal 271 end if 272 273 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 274 275 !Save initial values 276 if (itime==1)then 277 acell0(:)=acell(:) 278 rprimd0(:,:)=rprimd(:,:) 279 ucvol0=ucvol 280 end if 281 282 !DEBUG (UCVOL) 283 if(DEBUG)then 284 write(std_out,*) 'Volume of cell (ucvol):',ucvol 285 end if 286 287 !Get rid of mean force on whole unit cell, but only if no 288 !generalized constraints are in effect 289 ! residual_corrected(:,:)=residual(:,:) 290 ! if(ab_mover%nconeq==0)then 291 ! do ii=1,3 292 ! if (ii/=3.or.ab_mover%jellslab==0) then 293 ! favg=sum(residual_corrected(ii,:))/dble(ab_mover%natom) 294 ! residual_corrected(ii,:)=residual_corrected(ii,:)-favg 295 ! end if 296 ! end do 297 ! end if 298 299 !write(std_out,*) 'delocint 05' 300 !########################################################## 301 !### 05. Compute internals for first time 302 303 if (itime==1)then 304 call make_prim_internals(ab_mover%deloc,icenter,ab_mover%natom,& 305 & ab_mover%ntypat,rprimd,ab_mover%typat,xcart,ab_mover%znucl) 306 307 ABI_ALLOCATE(prim_int,(ab_mover%deloc%ninternal)) 308 309 if(DEBUG)then 310 write (message,'(a,i6)') 'Number of primitive internal coordinates (ninternal): ',ab_mover%deloc%ninternal 311 call wrtout(std_out, message,'COLL') 312 end if 313 314 if (allocated(u_matrix)) then 315 ABI_DEALLOCATE(u_matrix) 316 end if 317 ABI_ALLOCATE(u_matrix,(ab_mover%deloc%ninternal,ndeloc)) 318 319 call calc_prim_int(ab_mover%deloc,ab_mover%natom,rprimd,xcart,prim_int) 320 321 if(DEBUG)then 322 write (message,'(a)') 'Primitive internal coordinate values:' 323 call wrtout(std_out, message,'COLL') 324 write (message,'(a)') ' Bonds:' 325 call wrtout(std_out, message,'COLL') 326 do ii = 1, ab_mover%deloc%nbond 327 write (message,'(i6,E20.10)') ii, prim_int(ii) 328 call wrtout(std_out, message,'COLL') 329 end do 330 331 write (message,'(a)') ' Angles:' 332 call wrtout(std_out, message,'COLL') 333 do ii = ab_mover%deloc%nbond+1, ab_mover%deloc%nbond+ab_mover%deloc%nang 334 write (message,'(i6,2(E20.10,2x))') ii, prim_int(ii), prim_int(ii)/pi*180.0_dp 335 call wrtout(std_out, message,'COLL') 336 end do 337 338 write (message,'(a)') ' Dihedrals:' 339 call wrtout(std_out, message,'COLL') 340 do ii = ab_mover%deloc%nbond+ab_mover%deloc%nang+1, ab_mover%deloc%nbond+ab_mover%deloc%nang+ab_mover%deloc%ndihed 341 write (message,'(i6,2(E20.10,2x))') ii, prim_int(ii), prim_int(ii)/pi*180.0_dp 342 call wrtout(std_out, message,'COLL') 343 end do 344 345 write (message,'(a)') ' Cartesian auxiliary coordinates for constraints:' 346 call wrtout(std_out, message,'COLL') 347 do ii = ab_mover%deloc%nbond+ab_mover%deloc%nang+ab_mover%deloc%ndihed+1, ab_mover%deloc%ninternal 348 write (message,'(i6,E20.10)') ii, prim_int(ii) 349 call wrtout(std_out, message,'COLL') 350 end do 351 end if 352 353 ABI_DEALLOCATE(prim_int) 354 355 ! equal weight on all internal coordinates as a starting point. 356 u_matrix(:,:) = one / dble (ndeloc) 357 358 ! Zero the arrays before first use 359 deloc_force(:) = zero 360 361 end if 362 363 ABI_ALLOCATE(prim_int,(ab_mover%deloc%ninternal)) 364 365 !write(std_out,*) 'delocint 06' 366 !########################################################## 367 !### 06. Compute delocalized coordinates and forces 368 369 !xcart ---> deloc_int 370 371 !Convert positions to delocalized coordinates for next step 372 call xcart2deloc(ab_mover%deloc,ab_mover%natom,rprimd,xcart,& 373 & bt_inv_matrix,u_matrix,deloc_int,prim_int) 374 375 !fred ---> deloc_force 376 377 !Convert forces to delocalized coordinates for next step 378 call fred2fdeloc(bt_inv_matrix,deloc_force,residual,ab_mover%natom,gprimd) 379 380 !write(std_out,*) 'delocint 07' 381 !########################################################## 382 !### 07. Fill the vectors vin and vout 383 384 !DEBUG deloc_int and deloc_force before pack 385 if(DEBUG)then 386 write (std_out,*) 'Delocalized internals and forces (ndeloc):',ndeloc 387 write(std_out,*) 'deloc_int' 388 do ii=1,ndeloc,3 389 if (ii+2<=ndeloc)then 390 write(std_out,*) ii,deloc_int(ii:ii+2) 391 else 392 write(std_out,*) ii,deloc_int(ii:ndeloc) 393 end if 394 end do 395 write(std_out,*) 'deloc_force' 396 do ii=1,ndeloc,3 397 if (ii+2<=ndeloc)then 398 write(std_out,*) ii,deloc_force(ii:ii+2) 399 else 400 write(std_out,*) ii,deloc_force(ii:ndeloc) 401 end if 402 end do 403 end if 404 405 !DELOCINT 406 !Instead of fred_corrected we use deloc_force 407 !Instead of xred e use deloc_int 408 ! 409 !Initialize input vectors : first vin, then vout 410 !The values of vin from the previous iteration 411 !should be the same 412 call xfpack_x2vin(acell, acell0, ab_mover%natom-1, ndim,& 413 & ab_mover%nsym, ab_mover%optcell, rprimd, rprimd0,& 414 & ab_mover%symrel, ucvol, ucvol0, vin, deloc_int) 415 !end if 416 417 call xfpack_f2vout(deloc_force, ab_mover%natom-1, ndim,& 418 & ab_mover%optcell, ab_mover%strtarget, strten, ucvol, & 419 & vout) 420 421 !write(std_out,*) 'delocint 08' 422 !########################################################## 423 !### 08. Initialize or update the hessian matrix 424 425 !Initialise the Hessian matrix using gmet 426 if (itime==1)then 427 428 ! Initialise the Hessian matrix with ab_mover%userrc. 429 ! this has become unusable because it imposes ndim >= 3 natom 430 ! ident = 3x3 identity matrix 431 ! call hessinit(ab_mover, hessin, gmet, ndim, ucvol) 432 hessin = zero 433 do ii=1, ndim 434 hessin (ii,ii) = one 435 end do 436 437 ! ! Initialize inverse hessian with identity matrix 438 ! ! in cartesian coordinates, which makes use of metric tensor gmet 439 ! ! in reduced coordinates. 440 ! hessin(:,:)=zero 441 ! do ii=1,ab_mover%natom 442 ! do kk=1,3 443 ! do jj=1,3 444 ! ! Warning : implemented in reduced coordinates 445 ! if (ab_mover%iatfix(kk,ii)==0 .and.& 446 ! & ab_mover%iatfix(jj,ii)==0 )then 447 ! hessin(kk+3*(ii-1),jj+3*(ii-1))=gmet(kk,jj) 448 ! end if 449 ! end do 450 ! end do 451 ! end do 452 ! if(ab_mover%optcell/=0)then 453 ! ! These values might lead to too large changes in some cases 454 ! diag=ab_mover%strprecon*30.0_dp/ucvol 455 ! if(ab_mover%optcell==1) diag=diag/three 456 ! do ii=3*ab_mover%natom+1,ndim 457 ! hessin(ii,ii)=diag 458 ! end do 459 ! end if 460 461 if (ab_mover%restartxf/=0) then 462 463 call xfh_recover_deloc(ab_xfh,ab_mover,acell,acell0,cycl_main,& 464 & residual,hessin,ndim,rprimd,rprimd0,strten,ucvol,ucvol0,vin,vin_prev,& 465 & vout,vout_prev,xred,ab_mover%deloc,deloc_int,deloc_force,bt_inv_matrix,gprimd,prim_int,& 466 & u_matrix) 467 468 end if 469 470 end if 471 472 ABI_DEALLOCATE(prim_int) 473 474 if(itime>1)then 475 ! Update the hessian matrix, by taking into account the 476 ! current pair (x,f) and the previous one. 477 call hessupdt(hessin,ab_mover%iatfix,ab_mover%natom,ndim,vin,& 478 & vin_prev,vout,vout_prev) 479 480 end if 481 482 !DEBUG (vin,vout and hessin before prediction) 483 if(DEBUG)then 484 write(std_out,*) 'Vectors vin and vout and inverse of Hessian (hessin) [before]' 485 write(std_out,*) 'vin:' 486 do ii=1,ndim,3 487 if (ii+2<=ndim)then 488 write(std_out,*) ii,vin(ii:ii+2) 489 else 490 write(std_out,*) ii,vin(ii:ndim) 491 end if 492 end do 493 write(std_out,*) 'vout:' 494 do ii=1,ndim,3 495 if (ii+2<=ndim)then 496 write(std_out,*) ii,vout(ii:ii+2) 497 else 498 write(std_out,*) ii,vout(ii:ndim) 499 end if 500 end do 501 write(std_out,*) 'Inverse Hessian (hessin): ',ndim,'x',ndim 502 do kk=1,ndim 503 do jj=1,ndim,3 504 if (jj+2<=ndim)then 505 write(std_out,*) jj,hessin(jj:jj+2,kk) 506 else 507 write(std_out,*) jj,hessin(jj:ndim,kk) 508 end if 509 end do 510 end do 511 end if 512 513 !write(std_out,*) 'delocint 09' 514 !########################################################## 515 !### 09. Compute the next values 516 517 if(ionmov==10 .or. itime==1)then 518 519 ! Previous cartesian coordinates 520 vin_prev(:)=vin(:) 521 522 ! New atomic cartesian coordinates are obtained from vin, hessin 523 ! and vout 524 do ii=1,ndim 525 vin(:)=vin(:)-hessin(:,ii)*vout(ii) 526 end do 527 ! Previous atomic forces 528 vout_prev(:)=vout(:) 529 530 else 531 if(ionmov==11)then 532 ihist_prev = abihist_findIndex(hist,-1) 533 etotal_prev=hist%etot(ihist_prev) 534 ! Here the BFGS algorithm, modified to take into account the 535 ! energy 536 call brdene(etotal,etotal_prev,hessin,& 537 & ndim,vin,vin_prev,vout,vout_prev) 538 539 end if 540 541 ! DEBUG (vin,vout and hessin after prediction) 542 if(DEBUG)then 543 write(std_out,*) 'Vectors vin and vout [after prediction]' 544 write(std_out,*) 'vin:' 545 do ii=1,ndim,3 546 if (ii+2<=ndim)then 547 write(std_out,*) ii,vin(ii:ii+2) 548 else 549 write(std_out,*) ii,vin(ii:ndim) 550 end if 551 end do 552 write(std_out,*) 'vout:' 553 do ii=1,ndim,3 554 if (ii+2<=ndim)then 555 write(std_out,*) ii,vout(ii:ii+2) 556 else 557 write(std_out,*) ii,vout(ii:ndim) 558 end if 559 end do 560 end if 561 562 ! Implement fixing of atoms : put back old values for fixed 563 ! components 564 do kk=1,ab_mover%natom 565 do jj=1,3 566 ! Warning : implemented in reduced coordinates 567 if ( ab_mover%iatfix(jj,kk)==1) then 568 vin(jj+(kk-1)*3)=vin_prev(jj+(kk-1)*3) 569 end if 570 end do 571 end do 572 end if 573 574 !write(std_out,*) 'delocint 10' 575 !########################################################## 576 !### 10. Convert from delocalized to xcart and xred 577 578 !Transfer vin to deloc_int, acell and rprimd 579 call xfpack_vin2x(acell, acell0, ab_mover%natom-1, ndim,& 580 & ab_mover%nsym, ab_mover%optcell, rprimd, rprimd0,& 581 & ab_mover%symrel, ucvol, ucvol0,& 582 & vin, deloc_int) 583 584 if(DEBUG)then 585 write (std_out,*) 'Delocalized internals (deloc_int) [after prediction]:' 586 write(std_out,*) 'deloc_int:' 587 do ii=1,ndeloc,3 588 if (ii+2<=ndeloc)then 589 write(std_out,*) ii,deloc_int(ii:ii+2) 590 else 591 write(std_out,*) ii,deloc_int(ii:ndeloc) 592 end if 593 end do 594 write(std_out,*) 'BT Inverse Matrix:' 595 do ii=1,3*ab_mover%natom 596 write(std_out,*) bt_inv_matrix(:,ii) 597 end do 598 write (std_out,*) 'xcart (before deloc2xcart):' 599 do ii=1,ab_mover%natom 600 write (std_out,*) xcart(:,ii) 601 end do 602 end if 603 604 !this routine contains an iterative scheme to find xcart 605 !from the non-linear relations between deloc and xcart 606 !SIGNIFICANTLY DIFFERENT FROM xcart2deloc 607 call deloc2xcart(ab_mover%deloc,ab_mover%natom,rprimd,xcart,& 608 & deloc_int,bt_inv_matrix,u_matrix) 609 610 !Convert new xcart (cartesian) to xred (reduced coordinates) 611 call xcart2xred(ab_mover%natom,rprimd,xcart,xred) 612 613 614 !write(std_out,*) 'delocint 11' 615 !########################################################## 616 !### 11. Update the history with the prediction 617 618 !Increase indexes 619 hist%ihist = abihist_findIndex(hist,+1) 620 621 if(ab_mover%optcell/=0)then 622 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 623 end if 624 625 !Fill the history with the variables 626 !xred, acell, rprimd, vel 627 call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 628 ihist_prev = abihist_findIndex(hist,-1) 629 hist%vel(:,:,hist%ihist)=hist%vel(:,:,ihist_prev) 630 631 if(zDEBUG)then 632 write (std_out,*) 'residual:' 633 do kk=1,ab_mover%natom 634 write (std_out,*) residual(:,kk) 635 end do 636 write (std_out,*) 'strten:' 637 write (std_out,*) strten(1:3),ch10,strten(4:6) 638 write (std_out,*) 'etotal:' 639 write (std_out,*) etotal 640 end if 641 642 end subroutine pred_delocint