TABLE OF CONTENTS
- ABINIT/gwls_polarisability
- m_hamiltonian/epsilon_k
- m_hamiltonian/matrix_function_epsilon_k
- m_hamiltonian/Pk
- m_hamiltonian/set_dielectric_function_frequency
ABINIT/gwls_polarisability [ Modules ]
NAME
gwls_polarisability
FUNCTION
.
COPYRIGHT
Copyright (C) 2009-2018 ABINIT group (JLJ, BR, MC) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
PARENTS
CHILDREN
SOURCE
21 #if defined HAVE_CONFIG_H 22 #include "config.h" 23 #endif 24 25 #include "abi_common.h" 26 27 28 module gwls_polarisability 29 ! local modules 30 use m_gwls_utility 31 use gwls_wf 32 use gwls_valenceWavefunctions 33 use gwls_hamiltonian 34 use gwls_lineqsolver 35 36 ! abinit modules 37 use defs_basis 38 use m_errors 39 use m_profiling_abi 40 use m_bandfft_kpt 41 42 implicit none 43 save 44 private
m_hamiltonian/epsilon_k [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
epsilon_k
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_polarisability
CHILDREN
epsilon_k
SOURCE
679 subroutine epsilon_k(psi_out,psi_in,omega) 680 681 682 !This section has been created automatically by the script Abilint (TD). 683 !Do not modify the following lines by hand. 684 #undef ABI_FUNC 685 #define ABI_FUNC 'epsilon_k' 686 !End of the abilint section 687 688 implicit none 689 690 real(dp), intent(out) :: psi_out(2,npw_k) 691 real(dp), intent(in) :: psi_in(2,npw_k), omega(2) 692 693 ! ************************************************************************* 694 695 psi_out = psi_in 696 call sqrt_vc_k(psi_out) 697 call Pk(psi_out,omega) 698 call sqrt_vc_k(psi_out) 699 700 psi_out = psi_in - psi_out 701 end subroutine epsilon_k
m_hamiltonian/matrix_function_epsilon_k [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
matrix_function_epsilon_k
FUNCTION
.
INPUTS
OUTPUT
PARENTS
CHILDREN
epsilon_k
SOURCE
762 subroutine matrix_function_epsilon_k(vector_out,vector_in,Hsize) 763 !---------------------------------------------------------------------------------------------------- 764 ! This function is a simple wrapper around epsilon_k to be fed to the Lanczos 765 ! algorithm. 766 !---------------------------------------------------------------------------------------------------- 767 768 !This section has been created automatically by the script Abilint (TD). 769 !Do not modify the following lines by hand. 770 #undef ABI_FUNC 771 #define ABI_FUNC 'matrix_function_epsilon_k' 772 !End of the abilint section 773 774 implicit none 775 integer, intent(in) :: Hsize 776 complex(dpc), intent(out) :: vector_out(Hsize) 777 complex(dpc), intent(in) :: vector_in(Hsize) 778 779 780 ! local variables 781 real(dp) :: psik (2,npw_k) 782 real(dp) :: psik2(2,npw_k) 783 784 ! ************************************************************************* 785 786 ! convert from one format to the other 787 psik(1,:) = dble (vector_in(:)) 788 psik(2,:) = dimag(vector_in(:)) 789 790 ! act on vector 791 792 793 call epsilon_k(psik2 ,psik, matrix_function_omega) 794 795 ! convert back 796 vector_out = cmplx_1*psik2(1,:)+cmplx_i*psik2(2,:) 797 798 end subroutine matrix_function_epsilon_k
m_hamiltonian/Pk [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
Pk
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_polarisability
CHILDREN
epsilon_k
SOURCE
89 subroutine Pk(psi_inout,omega) 90 !=============================================================================== 91 ! 92 ! This routine applies the polarizability operator to an arbitrary state psi_inout. 93 ! 94 ! 95 ! A note about parallelism: 96 ! ------------------------- 97 ! The input/output is in "linear algebra" configuration, which is to say 98 ! that ALL processors have a fraction of the G-vectors for this state. Internally, 99 ! this routine will parallelise over bands and FFT, thus using the "FFT" 100 ! configuration of the data. For more information, see the lobpcgwf.F90, and 101 ! the article by F. Bottin et al. 102 !=============================================================================== 103 104 !This section has been created automatically by the script Abilint (TD). 105 !Do not modify the following lines by hand. 106 #undef ABI_FUNC 107 #define ABI_FUNC 'Pk' 108 use interfaces_14_hidewrite 109 use interfaces_18_timing 110 !End of the abilint section 111 112 implicit none 113 114 real(dp), intent(inout) :: psi_inout(2,npw_k) 115 real(dp), intent(in) :: omega(2) 116 117 logical :: omega_imaginary 118 119 integer :: v, mb, iblk 120 121 real(dp):: norm_omega 122 real(dp), allocatable :: psik(:,:) 123 real(dp), allocatable :: psik_alltoall(:,:), psik_wrk_alltoall(:,:) 124 real(dp), allocatable :: psik_in_alltoall(:,:), psik_tmp_alltoall(:,:) 125 126 real(dp), allocatable :: psik_ext(:,:), psik_ext_alltoall(:,:) 127 128 129 real(dp), allocatable :: psir(:,:,:,:), psir_ext(:,:,:,:) 130 131 integer :: cplex 132 integer :: recy_i 133 134 135 integer :: mpi_band_rank 136 137 real(dp) :: list_SQMR_frequencies(2) 138 real(dp) :: list_QMR_frequencies(2,2) 139 140 141 integer, save :: icounter = 0 142 real(dp), save :: total_time1, total_time2, total_time 143 144 145 integer :: num_op_v, i_op_v, case_op_v 146 real(dp):: factor_op_v 147 real(dp):: lbda 148 real(dp):: zz(2) 149 150 character(len=500) :: message 151 152 real(dp) :: tsec(2) 153 integer :: GWLS_TIMAB, OPTION_TIMAB 154 155 ! ************************************************************************* 156 157 GWLS_TIMAB = 1524 158 OPTION_TIMAB = 1 159 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 160 161 162 163 icounter = icounter + 1 164 call cpu_time(total_time1) 165 166 167 !======================================== 168 ! Allocate work arrays and define 169 ! important parameters 170 !======================================== 171 GWLS_TIMAB = 1525 172 OPTION_TIMAB = 1 173 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 174 175 176 cplex = 2 ! complex potential 177 178 mpi_band_rank = mpi_enreg%me_band 179 180 ABI_ALLOCATE(psik, (2,npw_kb)) 181 ABI_ALLOCATE(psik_alltoall, (2,npw_g)) 182 ABI_ALLOCATE(psik_wrk_alltoall, (2,npw_g)) 183 ABI_ALLOCATE(psik_tmp_alltoall, (2,npw_g)) 184 ABI_ALLOCATE(psik_in_alltoall, (2,npw_g)) 185 186 ABI_ALLOCATE(psir, (2,n4,n5,n6)) 187 ABI_ALLOCATE(psir_ext,(2,n4,n5,n6)) 188 189 190 OPTION_TIMAB = 2 191 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 192 193 194 !-------------------------------------------------------------------------------- 195 ! 196 ! The polarizability acting on a state | PSI > is given by 197 ! 198 ! Pk | PSI > = 2 sum_{v} | h_v^* phi_v > ( factor of 2 comes from sum on spin) 199 ! 200 ! | h_v > = Pc . OPERATOR_v . Pc | PSI^* phi_v > 201 ! 202 ! There are multiple cases to consider: 203 ! 204 ! CASE: 205 ! 1) |omega| = 0 206 ! OPERATOR_v = 1/[H-ev] 207 ! prefactor = -4 208 ! 209 ! 2) omega = i lbda (imaginary) 210 ! OPERATOR_v = [H-ev]/(lbda^2+[H-ev]^2) 211 ! prefactor = -4 212 ! 213 ! 3) omega = lbda (real) (USE SQMR) 214 ! OPERATOR_v = { 1/[H-ev+lbda]+ 1/[H-ev-lbda] } 215 ! prefactor = -2 216 ! 217 ! 4) omega = lbda (real) (USE QMR) 218 ! OPERATOR_v = { 1/[H-ev+lbda]+ 1/[H-ev-lbda] } 219 ! prefactor = -2 220 ! 221 ! It simplifies the code below to systematize the algorithm. 222 ! 223 !-------------------------------------------------------------------------------- 224 225 ! Check which part of omega is non-zero. Default is that omega is real. 226 if (abs(omega(2)) < 1.0d-12) then 227 omega_imaginary=.false. 228 norm_omega = omega(1) 229 230 elseif (abs(omega(1)) < 1.0d-12 .and. abs(omega(2)) > 1.0d-12) then 231 omega_imaginary = .true. 232 norm_omega = omega(2) 233 234 else 235 write(message,"(a,es16.8,3a)")& 236 "omega=",omega,",",ch10,& 237 "but either it's real or imaginary part need to be 0 for the polarisability routine to work." 238 MSG_ERROR(message) 239 end if 240 241 242 !----------------------------------------------------------------- 243 ! I) Prepare global values depending on the CASE being considered 244 !----------------------------------------------------------------- 245 246 if (norm_omega < 1.0D-12 ) then 247 case_op_v = 1 248 num_op_v = 1 249 factor_op_v = -4.0_dp 250 251 else if (omega_imaginary) then 252 case_op_v = 2 253 num_op_v = 1 254 factor_op_v = -4.0_dp 255 256 else if( .not. activate_inf_shift_poles) then 257 case_op_v = 3 258 num_op_v = 2 259 factor_op_v = -2.0_dp 260 261 else 262 case_op_v = 4 263 num_op_v = 2 264 factor_op_v = -2.0_dp 265 266 inf_shift_poles = dtset%zcut 267 268 269 write(message,*) " inf_shift_poles = ",inf_shift_poles 270 call wrtout(std_out,message,'COLL') 271 end if 272 273 274 write(message,10)" " 275 call wrtout(std_out,message,'COLL') 276 write(message,10) " Pk: applying the polarizability on states" 277 call wrtout(std_out,message,'COLL') 278 write(message,10) " ==============================================================" 279 call wrtout(std_out,message,'COLL') 280 write(message,12) " CASE : ",case_op_v 281 call wrtout(std_out,message,'COLL') 282 283 284 285 !----------------------------------------------------------------- 286 ! II) copy conjugate of initial wavefunction in local array, 287 ! and set inout array to zero. Each FFT row of processors 288 ! must have a copy of the initial wavefunction in FFT 289 ! configuration! 290 !----------------------------------------------------------------- 291 call cpu_time(time1) 292 293 GWLS_TIMAB = 1526 294 OPTION_TIMAB = 1 295 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 296 297 298 ABI_ALLOCATE(psik_ext,(2,npw_kb)) 299 ABI_ALLOCATE(psik_ext_alltoall,(2,npw_g)) 300 301 ! fill the array psik_ext with copies of the external state 302 do mb = 1, blocksize 303 psik_ext(:,(mb-1)*npw_k+1:mb*npw_k) = psi_inout(:,:) 304 end do 305 306 ! change configuration of the data, from LA to FFT 307 call wf_block_distribute(psik_ext, psik_ext_alltoall,1) ! LA -> FFT 308 ! Now every row of FFT processors has a copy of the external state. 309 310 ! Copy the external state to the real space format, appropriate for real space products to be 311 ! used later. 312 313 call g_to_r(psir_ext,psik_ext_alltoall) 314 psir_ext(2,:,:,:) = -psir_ext(2,:,:,:) 315 316 317 ! Don't need these arrays anymore... 318 ABI_DEALLOCATE(psik_ext) 319 ABI_DEALLOCATE(psik_ext_alltoall) 320 321 OPTION_TIMAB = 2 322 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 323 324 call cpu_time(time2) 325 time_fft = time_fft + time2-time1 326 counter_fft = counter_fft+1 327 328 ! set external state to zero, ready to start cumulating the answer 329 psi_inout = zero 330 331 !----------------------------------------------------------------- 332 ! III) Iterate on all valence bands 333 !----------------------------------------------------------------- 334 335 ! Loop on all blocks of eigenstates 336 do iblk = 1, nbdblock 337 338 339 ! What is the valence band index for this block and this row of FFT processors? It is not clear from the 340 ! code in lobpcgwf.F90; I'm going to *guess*. 341 342 v = (iblk-1)*blocksize + mpi_band_rank + 1 ! CAREFUL! This is a guess. Revisit this if code doesn't work as intended. 343 344 345 !Solving of Sternheiner equation 346 write(message,12) " band :", v 347 call wrtout(std_out,message,'COLL') 348 write(message,14) " eigenvalue (Ha) : ",eig(v) 349 call wrtout(std_out,message,'COLL') 350 write(message,14) " Re[omega] (Ha) : ",omega(1) 351 call wrtout(std_out,message,'COLL') 352 write(message,14) " Im[omega] (Ha) : ",omega(2) 353 call wrtout(std_out,message,'COLL') 354 355 !----------------------------------------------------------------- 356 ! IV) prepare some arrays, if they are needed 357 !----------------------------------------------------------------- 358 if (case_op_v == 3) then 359 360 list_SQMR_frequencies(1) = eig(v) - norm_omega 361 list_SQMR_frequencies(2) = eig(v) + norm_omega 362 363 else if (case_op_v == 4) then 364 list_QMR_frequencies(:,1) = (/eig(v)-norm_omega,-inf_shift_poles/) 365 list_QMR_frequencies(:,2) = (/eig(v)+norm_omega, inf_shift_poles/) 366 end if 367 368 369 !----------------------------------------------------------------- 370 ! V) Compute the real-space product of the input wavefunction 371 ! with the valence wavefunction. 372 !----------------------------------------------------------------- 373 ! Unfortunately, the input wavefunctions will be FFT-transformed k-> r nbandv times 374 ! because fourwf cannot multiply a real-space potential with a real-space wavefunction! 375 call cpu_time(time1) 376 377 GWLS_TIMAB = 1527 378 OPTION_TIMAB = 1 379 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 380 381 call gr_to_g(psik_in_alltoall,psir_ext,valence_wavefunctions_FFT(:,:,iblk)) 382 383 OPTION_TIMAB = 2 384 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 385 386 387 call cpu_time(time2) 388 time_fft = time_fft + time2-time1 389 counter_fft = counter_fft+1 390 391 !----------------------------------------------------------------- 392 ! VI) Project out to conduction space 393 !----------------------------------------------------------------- 394 call cpu_time(time1) 395 396 GWLS_TIMAB = 1528 397 OPTION_TIMAB = 1 398 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 399 400 !call pc_k(psik_in) 401 call pc_k_valence_kernel(psik_in_alltoall) 402 403 OPTION_TIMAB = 2 404 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 405 406 call cpu_time(time2) 407 time_proj = time_proj + time2-time1 408 counter_proj= counter_proj+1 409 410 411 !----------------------------------------------------------------- 412 ! VII) Loop on potential valence operators, 413 !----------------------------------------------------------------- 414 do i_op_v = 1, num_op_v 415 416 if (case_op_v == 1) then 417 418 ! frequency is zero, operator to apply is 1/[H-ev] 419 call cpu_time(time1) 420 GWLS_TIMAB = 1529 421 OPTION_TIMAB = 1 422 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 423 424 call sqmr(psik_in_alltoall,psik_tmp_alltoall,eig(v),1) 425 426 OPTION_TIMAB = 2 427 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 428 429 430 call cpu_time(time2) 431 time_sqmr = time_sqmr+ time2-time1 432 counter_sqmr= counter_sqmr+1 433 434 ! If we are constructing the $\hat \epsilon(i\omega = 0)$ matrix (and the Lanczos basis at the same time), 435 ! keep the Sternheimer solutions for use in the projected Sternheimer section (in LA configuration). 436 if(write_solution .and. ((iblk-1)*blocksize < nbandv)) then 437 call wf_block_distribute(psik, psik_tmp_alltoall, 2) ! FFT -> LA 438 do mb=1,blocksize 439 v = (iblk-1)*blocksize+mb 440 if(v <= nbandv) then 441 if(dtset%gwls_recycle == 1) then 442 Sternheimer_solutions_zero(:,:,index_solution,v) = psik(:,(mb-1)*npw_k+1:mb*npw_k) 443 end if 444 if(dtset%gwls_recycle == 2) then 445 recy_i = (index_solution-1)*nbandv + v 446 !BUG : On petrus, NAG 5.3.1 + OpenMPI 1.6.2 cause read(...,rec=i) to read the data written by write(...,rec=i+1). 447 !Workaround compatible only with nag : write(recy_unit,rec=recy_i+1). 448 write(recy_unit,rec=recy_i) psik(:,(mb-1)*npw_k+1:mb*npw_k) 449 end if 450 end if 451 end do 452 end if 453 454 else if (case_op_v == 2) then 455 456 ! frequency purely imaginary, operator to apply is 457 ! [H-ev]/(lbda^2+[H-ev]^2) 458 call cpu_time(time1) 459 460 GWLS_TIMAB = 1533 461 OPTION_TIMAB = 1 462 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 463 464 psik_wrk_alltoall = psik_in_alltoall 465 call Hpsik(psik_wrk_alltoall,eig(v)) 466 467 OPTION_TIMAB = 2 468 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 469 470 471 call cpu_time(time2) 472 time_H = time_H + time2-time1 473 counter_H = counter_H+1 474 475 call cpu_time(time1) 476 477 GWLS_TIMAB = 1530 478 OPTION_TIMAB = 1 479 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 480 481 call sqmr(psik_wrk_alltoall,psik_tmp_alltoall,eig(v),0,norm_omega,omega_imaginary) 482 483 OPTION_TIMAB = 2 484 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 485 486 487 call cpu_time(time2) 488 time_sqmr = time_sqmr+ time2-time1 489 counter_sqmr= counter_sqmr+1 490 491 else if (case_op_v == 3) then 492 493 ! frequency purely real, operator to apply is 494 ! 1/[H-ev +/- lbda] 495 ! TREATED WITH SQMR 496 497 498 lbda = list_SQMR_frequencies(i_op_v) 499 call cpu_time(time1) 500 501 GWLS_TIMAB = 1531 502 OPTION_TIMAB = 1 503 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 504 505 call sqmr(psik_in_alltoall,psik_tmp_alltoall,lbda,1) 506 507 OPTION_TIMAB = 2 508 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 509 510 511 512 call cpu_time(time2) 513 time_sqmr = time_sqmr+ time2-time1 514 counter_sqmr= counter_sqmr+1 515 516 else if (case_op_v == 4) then 517 ! frequency purely real, operator to apply is 518 ! 1/[H-ev +/- lbda] 519 ! TREATED WITH QMR 520 521 zz(:) = list_QMR_frequencies(:,i_op_v) 522 523 call cpu_time(time1) 524 GWLS_TIMAB = 1532 525 OPTION_TIMAB = 1 526 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 527 528 call qmr(psik_in_alltoall,psik_tmp_alltoall,zz) !,0) 529 530 OPTION_TIMAB = 2 531 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 532 533 534 call cpu_time(time2) 535 time_sqmr = time_sqmr+ time2-time1 536 counter_sqmr= counter_sqmr+1 537 538 end if 539 !----------------------------------------------------------------- 540 ! VIII) Project on conduction states 541 !----------------------------------------------------------------- 542 call cpu_time(time1) 543 GWLS_TIMAB = 1528 544 OPTION_TIMAB = 1 545 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 546 547 call pc_k_valence_kernel(psik_tmp_alltoall) 548 549 OPTION_TIMAB = 2 550 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 551 552 553 call cpu_time(time2) 554 time_proj = time_proj + time2-time1 555 counter_proj= counter_proj+1 556 557 !----------------------------------------------------------------- 558 ! IX) Conjugate result, and express in denpot format 559 !----------------------------------------------------------------- 560 561 call cpu_time(time1) 562 563 GWLS_TIMAB = 1526 564 OPTION_TIMAB = 1 565 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 566 567 568 call g_to_r(psir,psik_tmp_alltoall) 569 psir(2,:,:,:) = -psir(2,:,:,:) 570 571 572 OPTION_TIMAB = 2 573 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 574 575 call cpu_time(time2) 576 time_fft = time_fft + time2-time1 577 counter_fft = counter_fft+1 578 579 !----------------------------------------------------------------- 580 ! X) Multiply by valence state in real space 581 !----------------------------------------------------------------- 582 call cpu_time(time1) 583 GWLS_TIMAB = 1527 584 OPTION_TIMAB = 1 585 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 586 587 call gr_to_g(psik_alltoall,psir,valence_wavefunctions_FFT(:,:,iblk)) 588 589 590 OPTION_TIMAB = 2 591 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 592 593 594 call cpu_time(time2) 595 time_fft = time_fft + time2-time1 596 counter_fft = counter_fft+1 597 598 !----------------------------------------------------------------- 599 ! XI) Return to LA configuration, and cumulate the sum 600 !----------------------------------------------------------------- 601 call wf_block_distribute(psik, psik_alltoall,2) ! FFT -> LA 602 603 do mb = 1, blocksize 604 605 v = (iblk-1)*blocksize + mb 606 607 if (v <= nbandv) then 608 ! only add contributions from valence 609 psi_inout(:,:) = psi_inout(:,:) + psik(:,(mb-1)*npw_k+1:mb*npw_k) 610 end if 611 612 end do 613 614 615 end do ! i_op_v 616 617 end do ! iblk 618 619 620 !----------------------------------------------------------------- 621 ! XII) account for prefactor 622 !----------------------------------------------------------------- 623 psi_inout = factor_op_v*psi_inout 624 625 626 GWLS_TIMAB = 1525 627 OPTION_TIMAB = 1 628 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 629 630 631 ABI_DEALLOCATE(psik) 632 ABI_DEALLOCATE(psik_alltoall) 633 ABI_DEALLOCATE(psik_wrk_alltoall) 634 ABI_DEALLOCATE(psik_tmp_alltoall) 635 ABI_DEALLOCATE(psik_in_alltoall) 636 637 ABI_DEALLOCATE(psir) 638 ABI_DEALLOCATE(psir_ext) 639 640 OPTION_TIMAB = 2 641 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 642 643 644 call cpu_time(total_time2) 645 646 total_time = total_time + total_time2 - total_time1 647 648 649 GWLS_TIMAB = 1524 650 OPTION_TIMAB = 2 651 call timab(GWLS_TIMAB,OPTION_TIMAB,tsec) 652 653 10 format(A) 654 12 format(A,I5) 655 14 format(A,F24.16) 656 657 end subroutine Pk
m_hamiltonian/set_dielectric_function_frequency [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
set_dielectric_function_frequency
FUNCTION
.
INPUTS
OUTPUT
PARENTS
gwls_ComputeCorrelationEnergy,gwls_ComputePoles,gwls_GenerateEpsilon
CHILDREN
epsilon_k
SOURCE
723 subroutine set_dielectric_function_frequency(omega) 724 !---------------------------------------------------------------------------------------------------- 725 ! This routine sets the value of the module's frequency. 726 !---------------------------------------------------------------------------------------------------- 727 728 !This section has been created automatically by the script Abilint (TD). 729 !Do not modify the following lines by hand. 730 #undef ABI_FUNC 731 #define ABI_FUNC 'set_dielectric_function_frequency' 732 !End of the abilint section 733 734 implicit none 735 real(dp), intent(in) :: omega(2) 736 737 ! ************************************************************************* 738 739 matrix_function_omega(:) = omega(:) 740 741 end subroutine set_dielectric_function_frequency