Line data Source code
1 : !!****m* ABINIT/m_raman
2 : !! NAME
3 : !! m_raman
4 : !!
5 : !! FUNCTION
6 : !! Raman susceptibilities of zone-center phonons and electroo tensor.
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 1999-2026 ABINIT group (MVeithen)
10 : !! This file is distributed under the terms of the
11 : !! GNU General Public License, see ~abinit/COPYING
12 : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : !!
14 : !! SOURCE
15 :
16 : #if defined HAVE_CONFIG_H
17 : #include "config.h"
18 : #endif
19 :
20 : #include "abi_common.h"
21 :
22 : module m_raman
23 :
24 : use defs_basis
25 : use m_errors
26 : use m_abicore
27 : use m_nctk
28 : use netcdf
29 :
30 : use m_fstrings, only : sjoin
31 : use m_matrix, only : matr3inv
32 :
33 : implicit none
34 :
35 : private
36 : !!***
37 :
38 : public :: ramansus ! Raman susceptibilities of zone-center phonons.
39 : public :: electrooptic ! Electrooptic tensor and the raman tensors of zone-center phonons.
40 : public :: defwrite_nonana_raman_terms ! Write raman susceptiblities for q-->0
41 : public :: defwrite_raman_terms ! Write raman susceptiblities and frequencies for q=0
42 : !!***
43 :
44 : contains
45 : !!***
46 :
47 : !!****f* m_raman/ramansus
48 : !!
49 : !! NAME
50 : !! ramansus
51 : !!
52 : !! FUNCTION
53 : !! Compute the raman susceptibilities of zone-center phonons
54 : !!
55 : !! INPUTS
56 : !! d2cart = second order derivatives of the energy wrt all perturbations
57 : !! dchide(3,3,3) = non-linear optical coefficients from dtchi
58 : !! dchidt(natom,3,3,3) = first-order change of the electronic dielectric
59 : !! tensor induced by an individual atomic displacement
60 : !! displ = phonon mode atomic displacements
61 : !! mpert = maximum number of perturbations
62 : !! natom = number of atoms
63 : !! phfrq = phonon frequencies
64 : !! qphnrm=(described below)
65 : !! ucvol = unit cell volume
66 : !!
67 : !! OUTPUT
68 : !! qphon(3)= to be divided by qphnrm, give the phonon wavevector;
69 : !! if qphnrm==0.0_dp, then the wavevector is zero (Gamma point)
70 : !! and qphon gives the direction of
71 : !! the induced electric field; in the latter case, if qphon is
72 : !! zero, no non-analytical contribution is included.
73 : !! rsus
74 : !!
75 : !! SOURCE
76 :
77 6 : subroutine ramansus(d2cart,dchide,dchidt,displ,mpert,natom,phfrq,qphon,qphnrm,rsus,ucvol)
78 :
79 : !Arguments -----------------------------------
80 : !scalars
81 : integer,intent(in) :: mpert,natom
82 : real(dp),intent(in) :: qphnrm,ucvol
83 : !arrays
84 : real(dp),intent(in) :: d2cart(2,3,mpert,3,mpert),dchide(3,3,3)
85 : real(dp),intent(in) :: dchidt(natom,3,3,3),displ(2,3*natom,3*natom)
86 : real(dp),intent(in) :: phfrq(3*natom)
87 : real(dp),intent(inout) :: qphon(3)
88 : real(dp),intent(out) :: rsus(3*natom,3,3)
89 :
90 : !Local variables-------------------------------
91 : !scalars
92 : integer :: analyt,i1,i1dir,i1pert,i2dir,iatom,idir,imode
93 : real(dp) :: epsq,fac,g0,g1,g2,qphon2
94 : logical :: t_degenerate,iwrite
95 : character(len=500) :: message
96 : !arrays
97 : real(dp) :: dijk_q(3,3)
98 6 : real(dp),allocatable :: zeff(:,:)
99 6 : character(len=1),allocatable :: metacharacter(:)
100 :
101 : ! *********************************************************************
102 :
103 6 : iwrite = ab_out > 0
104 :
105 18 : ABI_MALLOC(zeff,(3,natom))
106 :
107 402 : rsus(:,:,:) = zero
108 6 : epsq = zero
109 54 : zeff(:,:) = zero
110 : dijk_q(:,:) = zero
111 :
112 : !Determine the analyticity of the matrix.
113 6 : analyt=1
114 6 : if(abs(qphnrm)<tol8)analyt=0
115 6 : if(abs(qphon(1))<tol8.and.abs(qphon(2))<tol8.and.abs(qphon(3))<tol8) analyt=1
116 :
117 : !In the case the non-analyticity is required :
118 4 : if(analyt == 0) then
119 :
120 : ! Normalize the limiting direction
121 4 : qphon2=qphon(1)**2+qphon(2)**2+qphon(3)**2
122 4 : qphon(1)=qphon(1)/sqrt(qphon2)
123 4 : qphon(2)=qphon(2)/sqrt(qphon2)
124 4 : qphon(3)=qphon(3)/sqrt(qphon2)
125 :
126 : ! Get the dielectric constant for the limiting direction
127 : epsq= 0._dp
128 16 : do i1dir=1,3
129 52 : do i2dir=1,3
130 48 : epsq=epsq+qphon(i1dir)*qphon(i2dir)*d2cart(1,i1dir,natom+2,i2dir,natom+2)
131 : end do
132 : end do
133 :
134 : ! Check if epsq > 0
135 4 : if (epsq < tol8) then
136 0 : write(message,'(a,es14.6)')' The value of epsq must be > 0 while it is found to be',epsq
137 0 : ABI_BUG(message)
138 : end if
139 :
140 : ! Get the effective charges for the limiting direction
141 16 : do i1dir=1,3
142 40 : do i1pert=1,natom
143 24 : zeff(i1dir,i1pert)=zero
144 108 : do i2dir=1,3
145 : zeff(i1dir,i1pert)=zeff(i1dir,i1pert)+qphon(i2dir)*&
146 96 : & d2cart(1,i1dir,i1pert,i2dir,natom+2)
147 : end do
148 : end do
149 : end do
150 :
151 : ! Get the NLO tensor for the limiting direction !$\sum_{k} d_{ijk} \cdot q_k$
152 :
153 4 : dijk_q(:,:) = zero
154 16 : do i1dir = 1, 3
155 52 : do i2dir = 1, 3
156 156 : do idir = 1, 3
157 144 : dijk_q(i1dir,i2dir) = dijk_q(i1dir,i2dir) + dchide(i1dir,i2dir,idir)*qphon(idir)
158 : end do
159 : end do
160 : end do
161 :
162 4 : fac = 16._dp*pi/(ucvol*epsq)
163 28 : do imode = 1, 3*natom
164 76 : do iatom = 1, natom
165 216 : do idir = 1, 3
166 144 : i1=idir + (iatom - 1)*3
167 : rsus(imode,:,:) = rsus(imode,:,:) + &
168 : & (dchidt(iatom,idir,:,:) - fac*zeff(idir,iatom)*dijk_q(:,:))* &
169 1920 : & displ(1,i1,imode)
170 : end do ! disp
171 : end do ! iatom
172 : end do ! imode
173 268 : rsus(:,:,:) = rsus(:,:,:)*sqrt(ucvol)
174 :
175 : else
176 14 : do imode = 1, 3*natom
177 38 : do iatom = 1, natom
178 108 : do idir = 1, 3
179 72 : i1=idir + (iatom - 1)*3
180 960 : rsus(imode,:,:) = rsus(imode,:,:) + dchidt(iatom,idir,:,:)*displ(1,i1,imode)
181 : end do ! disp
182 : end do ! iatom
183 : end do ! imode
184 134 : rsus(:,:,:) = rsus(:,:,:)*sqrt(ucvol)
185 : end if ! analyt == 0
186 :
187 6 : if (analyt == 0) then
188 4 : if (iwrite) then
189 4 : write(ab_out,*) ch10
190 : write(ab_out, '(a,/,a,3f9.5)' )&
191 4 : & ' Raman susceptibility of zone-center phonons, with non-analyticity in the',&
192 20 : & ' direction (cartesian coordinates)',qphon(1:3)+tol10
193 : write(ab_out,'(a)')&
194 4 : & ' -----------------------------------------------------------------------'
195 4 : write(ab_out,*) ch10
196 : end if
197 :
198 : else
199 2 : if (iwrite) then
200 2 : write(ab_out,*) ch10
201 2 : write(ab_out,*)' Raman susceptibilities of transverse zone-center phonon modes'
202 2 : write(ab_out,*)' -------------------------------------------------------------'
203 2 : write(ab_out,*) ch10
204 : end if
205 : end if
206 :
207 : !Examine the degeneracy of each mode. The portability of the echo of the Raman susceptibility
208 : !for each degenerate mode is very hard to guarantee. On the contrary,
209 : !the scalar reductions of these quantities are OK.
210 24 : ABI_MALLOC(metacharacter,(3*natom))
211 42 : do imode=1,3*natom
212 : ! The degenerate modes are not portable
213 36 : t_degenerate=.false.
214 36 : if(imode>1)then
215 30 : if(phfrq(imode)-phfrq(imode-1)<tol6)t_degenerate=.true.
216 : end if
217 36 : if(imode<3*natom)then
218 30 : if(phfrq(imode+1)-phfrq(imode)<tol6)t_degenerate=.true.
219 : end if
220 36 : metacharacter(imode)=';'
221 42 : if(t_degenerate)metacharacter(imode)='-'
222 : end do
223 :
224 42 : do imode = 1, 3*natom
225 36 : if (iwrite) then
226 36 : write(ab_out,'(a4,i3,2x,a2,f7.2,a6)')' Mode',imode,' (',phfrq(imode)*Ha_cmm1,' cm-1)'
227 144 : do idir = 1,3
228 144 : write(ab_out,'(a,4x,3(f16.9,2x))')metacharacter(imode),rsus(imode,idir,:)
229 : end do
230 : end if
231 : ! See R. Caracas and X. Gonze, Thermodynamic Properties of Solids : experiment and modeling, Wiley-VCH,
232 : ! Ed. S. Chaplot and R. Mittal and N. Choudhury , chap. 8, pp 291-312.
233 36 : g0=(rsus(imode,1,1)+rsus(imode,2,2)+rsus(imode,3,3))**2*third
234 : g1=((rsus(imode,1,2)-rsus(imode,2,1))**2+&
235 : & (rsus(imode,1,3)-rsus(imode,3,1))**2+&
236 36 : & (rsus(imode,2,3)-rsus(imode,3,2))**2)*half
237 : g2=((rsus(imode,1,2)+rsus(imode,2,1))**2+&
238 : & (rsus(imode,1,3)+rsus(imode,3,1))**2+&
239 : & (rsus(imode,2,3)+rsus(imode,3,2))**2)*half +&
240 : & ((rsus(imode,1,1)-rsus(imode,2,2))**2+&
241 : & (rsus(imode,2,2)-rsus(imode,3,3))**2+&
242 36 : & (rsus(imode,3,3)-rsus(imode,1,1))**2)*third
243 42 : if (iwrite) then
244 36 : write(ab_out,'(3(a,f16.9))')' Spherical averages : G0=',g0,' G1=',g1,' G2=',g2
245 36 : write(ab_out,*)
246 : end if
247 : end do
248 :
249 6 : ABI_FREE(metacharacter)
250 6 : ABI_FREE(zeff)
251 :
252 6 : end subroutine ramansus
253 : !!***
254 :
255 : !!****f* m_raman/electrooptic
256 : !! NAME
257 : !! electrooptic
258 : !!
259 : !! FUNCTION
260 : !! Compute the electrooptic tensor and the raman tensors of zone-center phonons
261 : !!
262 : !! INPUTS
263 : !! dchide(3,3,3) = non-linear optical coefficients
264 : !! dieflag= dielectric tensor flag. 0=> no dielectric tensor,
265 : !! 1=> frequency-dependent dielectric tensor,
266 : !! 2=> only the electronic dielectric tensor.
267 : !! epsinf=electronic dielectric tensor
268 : !! fact_oscstr(2,3,3*natom)=factors of the oscillator strengths for the different eigenmodes,
269 : !! for different direction of the electric field
270 : !! natom=number of atoms in unit cell
271 : !! phfrq(3*natom)=phonon frequencies (square root of the dynamical
272 : !! matrix eigenvalues, except if these are negative, and in this
273 : !! case, give minus the square root of the absolute value
274 : !! of the matrix eigenvalues). Hartree units.
275 : !! prtmbm= if equal to 1 write out the mode by mode decomposition of the EO tensor
276 : !! rsus = Raman susceptibilities
277 : !! ucvol=unit cell volume
278 : !!
279 : !! OUTPUT
280 : !! (to be completed ?)
281 : !!
282 : !! NOTES
283 : !! 1. The phonon frequencies phfrq should correspond to the
284 : !! wavevector at Gamma, without any non-analyticities.
285 : !! 2. Should clean for no imaginary part ...
286 : !! This routine should be used only by one processor.
287 : !! 3. frdiel(3,3,nfreq)= frequency-dependent dielectric tensor
288 : !! mode effective charges for the different eigenmodes,
289 : !! for different direction of the electric field
290 : !!
291 : !! SOURCE
292 :
293 2 : subroutine electrooptic(dchide,dieflag,epsinf,fact_oscstr,natom,phfrq,prtmbm,rsus,ucvol)
294 :
295 : !Arguments -------------------------------
296 : !scalars
297 : integer,intent(in) :: dieflag,natom,prtmbm
298 : real(dp),intent(in) :: ucvol
299 : !arrays
300 : real(dp),intent(in) :: dchide(3,3,3),epsinf(3,3),fact_oscstr(2,3,3*natom)
301 : real(dp),intent(in) :: phfrq(3*natom),rsus(3*natom,3,3)
302 :
303 : !Local variables -------------------------
304 : !scalars
305 : integer :: flag,i1,i2,ii,imode,jj,kk
306 : real(dp) :: dtm,fac
307 : logical :: iwrite
308 : character(len=500) :: message
309 : !arrays
310 : integer :: voigtindex(6,2)
311 : real(dp) :: eta(3,3),rvoigt(6,3),work(3,3,3)
312 2 : real(dp),allocatable :: rijk(:,:,:,:),rijk_tot(:,:,:)
313 :
314 : ! *********************************************************************
315 :
316 : !rijk(1:3*natom,:,:,:) = mode by mode decomposition of the electrooptic tensor
317 : !rijk(3*natom+1,:,:,:) = electronic contribution
318 2 : iwrite = ab_out > 0
319 :
320 2 : voigtindex(1,1) = 1 ; voigtindex(1,2) = 1
321 2 : voigtindex(2,1) = 2 ; voigtindex(2,2) = 2
322 2 : voigtindex(3,1) = 3 ; voigtindex(3,2) = 3
323 2 : voigtindex(4,1) = 2 ; voigtindex(4,2) = 3
324 2 : voigtindex(5,1) = 1 ; voigtindex(5,2) = 3
325 2 : voigtindex(6,1) = 1 ; voigtindex(6,2) = 2
326 :
327 12 : ABI_MALLOC(rijk,(3*natom+1,3,3,3))
328 2 : ABI_MALLOC(rijk_tot,(3,3,3))
329 458 : rijk(:,:,:,:) = 0._dp
330 80 : rijk_tot(:,:,:) = 0._dp
331 :
332 :
333 : !In case there is no mode with truly negative frequency
334 : !and the electronic dielectric tensor is available
335 : !compute the electro-optic tensor
336 :
337 2 : flag = 1
338 :
339 2 : if (abs(phfrq(1)) > abs(phfrq(4))) then
340 0 : flag = 0
341 : write(message,'(6a)')&
342 0 : & 'The lowest mode appears to be a "true" negative mode,',ch10,&
343 0 : & 'and not an acoustic mode. This precludes the computation',ch10,&
344 0 : & 'of the EO tensor.',ch10
345 0 : ABI_WARNING(message)
346 : end if
347 :
348 : dtm = epsinf(1,1)*epsinf(2,2)*epsinf(3,3) + &
349 : & epsinf(1,2)*epsinf(2,3)*epsinf(3,1) + &
350 : & epsinf(1,3)*epsinf(2,1)*epsinf(3,2) - &
351 : & epsinf(3,1)*epsinf(2,2)*epsinf(1,3) - &
352 : & epsinf(3,2)*epsinf(2,3)*epsinf(1,1) - &
353 2 : & epsinf(3,3)*epsinf(2,1)*epsinf(1,2)
354 :
355 2 : if (abs(dtm) < tol6) then
356 0 : flag = 0
357 : write(message,'(a,a,a,a,a,a,a,a)')&
358 0 : & 'The determinant of the electronic dielectric tensor is zero.',ch10,&
359 0 : & 'This preludes the computation fo the EO tensor since',ch10,&
360 0 : & 'this quantity requires the inverse of epsilon.',ch10,&
361 0 : & 'Action : check you database and the value of dieflag in the input file.',ch10
362 0 : ABI_WARNING(message)
363 : end if
364 :
365 : !dieflag is required to be one since the EO tensor
366 : !requires the oscillator strengths
367 :
368 2 : if ((flag == 1).and.(dieflag==1)) then
369 :
370 : ! Factor to convert atomic units to MKS units
371 :
372 2 : fac = -16._dp*pi*pi*eps0*(Bohr_Ang**2)*1.0d-8/(e_Cb*sqrt(ucvol))
373 :
374 : ! Compute inverse of dielectric tensor
375 : ! needed to convert the nonlinear optical susceptibility tensor
376 : ! to the electrooptic tensor
377 :
378 2 : call matr3inv(epsinf,eta)
379 :
380 2 : if (iwrite) then
381 2 : write(ab_out,*)ch10
382 2 : write(ab_out,*)'Output of the EO tensor (pm/V) in Voigt notations'
383 2 : write(ab_out,*)'================================================='
384 2 : write(ab_out,*)
385 2 : if (prtmbm == 1) then
386 2 : write(ab_out,*)'Mode by mode decomposition'
387 2 : write(ab_out,*)
388 : end if
389 : end if
390 :
391 : ! Compute the ionic contribution to the EO tensor
392 :
393 8 : do imode = 4, 3*natom
394 :
395 6 : if (prtmbm == 1 .and. iwrite) then
396 6 : write(ab_out,*)
397 6 : write(ab_out,'(a4,i3,2x,a2,f7.2,a6)')'Mode',imode,' (',phfrq(imode)*Ha_cmm1,' cm-1)'
398 : end if
399 :
400 24 : do ii = 1, 3
401 78 : do jj = 1, 3
402 234 : do kk = 1, 3
403 216 : rijk(imode,ii,jj,kk) = rsus(imode,ii,jj)*fact_oscstr(1,kk,imode)/(phfrq(imode)**2)
404 : end do
405 : end do
406 : end do
407 :
408 6 : work(:,:,:) = 0._dp
409 24 : do ii = 1,3
410 78 : do jj = 1, 3
411 234 : do kk = 1, 3
412 :
413 648 : do i1 = 1, 3
414 2106 : do i2 = 1, 3
415 1944 : work(ii,jj,kk) = work(ii,jj,kk) + eta(ii,i1)*rijk(imode,i1,i2,kk)*eta(i2,jj)
416 : end do ! i2
417 : end do ! i1
418 :
419 162 : rijk(imode,ii,jj,kk) = fac*work(ii,jj,kk)
420 216 : rijk_tot(ii,jj,kk) = rijk_tot(ii,jj,kk) + rijk(imode,ii,jj,kk)
421 : end do
422 :
423 : end do
424 : end do
425 :
426 8 : if (prtmbm == 1) then
427 6 : rvoigt(:,:) = 0._dp
428 42 : do i1 = 1, 6
429 36 : ii = voigtindex(i1,1)
430 36 : jj = voigtindex(i1,2)
431 144 : do kk = 1, 3
432 144 : rvoigt(i1,kk) = (rijk(imode,ii,jj,kk) + rijk(imode,jj,ii,kk))/2._dp
433 : end do
434 42 : if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:)
435 : end do
436 : end if
437 :
438 : end do ! imode
439 :
440 : ! Compute the electronic contribution to the EO tensor
441 :
442 2 : if (prtmbm == 1 .and. iwrite) then
443 2 : write(ab_out,*)
444 2 : write(ab_out,*)'Electronic contribution to the EO tensor'
445 : end if
446 :
447 8 : fac = 16*(pi**2)*(Bohr_Ang**2)*1.0d-8*eps0/e_Cb
448 :
449 8 : do ii = 1,3
450 26 : do jj = 1, 3
451 78 : do kk = 1, 3
452 :
453 216 : do i1 = 1, 3
454 702 : do i2 = 1, 3
455 : rijk(3*natom+1,ii,jj,kk) = rijk(3*natom+1,ii,jj,kk) + &
456 648 : & eta(ii,i1)*dchide(i1,i2,kk)*eta(i2,jj)
457 : end do ! i2
458 : end do ! i1
459 :
460 54 : rijk(3*natom+1,ii,jj,kk) = -4._dp*rijk(3*natom+1,ii,jj,kk)*fac
461 72 : rijk_tot(ii,jj,kk) = rijk_tot(ii,jj,kk) + rijk(3*natom+1,ii,jj,kk)
462 :
463 : end do
464 :
465 : end do
466 : end do
467 :
468 2 : if (prtmbm == 1) then
469 2 : rvoigt(:,:) = 0._dp
470 14 : do i1 = 1, 6
471 12 : ii = voigtindex(i1,1)
472 12 : jj = voigtindex(i1,2)
473 48 : do kk = 1, 3
474 48 : rvoigt(i1,kk) = (rijk(3*natom+1,ii,jj,kk) + rijk(3*natom+1,jj,ii,kk))/2._dp
475 : end do
476 14 : if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:)
477 : end do
478 2 : if (iwrite) write(ab_out,*)ch10
479 : end if
480 :
481 2 : if (iwrite) write(ab_out,*)'Total EO tensor (pm/V) in Voigt notations'
482 2 : rvoigt(:,:) = 0._dp
483 14 : do i1 = 1, 6
484 12 : ii = voigtindex(i1,1)
485 12 : jj = voigtindex(i1,2)
486 48 : do kk = 1, 3
487 48 : rvoigt(i1,kk) = (rijk_tot(ii,jj,kk) + rijk_tot(jj,ii,kk))/2._dp
488 : end do
489 14 : if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:)
490 : end do
491 :
492 : end if ! flag
493 :
494 2 : ABI_FREE(rijk)
495 2 : ABI_FREE(rijk_tot)
496 :
497 2 : end subroutine electrooptic
498 : !!***
499 :
500 : !!****f* m_raman/defwrite_nonana_raman_terms
501 : !! NAME
502 : !! defwrite_nonana_raman_terms
503 : !!
504 : !! FUNCTION
505 : !! Write the Raman susceptiblities for q-->0 along different directions in the netcdf file.
506 : !!
507 : !! INPUTS
508 : !! ncid=netcdf file id.
509 : !! iq_dir=Index of the q-point to be written to file.
510 : !! ndirs=Number of qpoints.
511 : !! rsus(3*natom,3,3)=List of Raman susceptibilities along the direction corresponding to iq_dir.
512 : !! natom=Number of atoms
513 : !!
514 : !! OUTPUT
515 : !! Only writing.
516 : !!
517 : !! SOURCE
518 :
519 6 : subroutine defwrite_nonana_raman_terms(ncid, iq_dir, ndirs, natom, rsus, mode)
520 :
521 : !Arguments ------------------------------------
522 : !scalars
523 : integer,intent(in) :: ncid,natom,iq_dir,ndirs
524 : character(len=*),intent(in) :: mode
525 : !arrays
526 : real(dp),intent(in) :: rsus(3*natom,3,3)
527 :
528 : !Local variables-------------------------------
529 : !scalars
530 : integer :: ncerr, raman_sus_varid
531 : ! *************************************************************************
532 :
533 : ! Fake use of ndirs, to keep it as argument. This should be removed when ndirs will be used.
534 : if(.false.) ncerr=ndirs
535 :
536 2 : select case (mode)
537 : case ("define")
538 2 : NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.))
539 : ncerr = nctk_def_arrays(ncid, [ nctkarr_t("non_analytical_raman_sus", "dp", &
540 4 : "number_of_non_analytical_directions,number_of_phonon_modes,number_of_cartesian_directions,number_of_cartesian_directions")])
541 2 : NCF_CHECK(ncerr)
542 :
543 2 : NCF_CHECK(nctk_set_datamode(ncid))
544 :
545 : case ("write")
546 4 : NCF_CHECK(nf90_inq_varid(ncid, "non_analytical_raman_sus", raman_sus_varid))
547 36 : ncerr = nf90_put_var(ncid,raman_sus_varid,rsus, start=[iq_dir,1,1,1], count=[1,3*natom,3,3])
548 4 : NCF_CHECK(ncerr)
549 :
550 : case default
551 6 : ABI_ERROR(sjoin("Wrong value for mode", mode))
552 : end select
553 :
554 6 : end subroutine defwrite_nonana_raman_terms
555 : !!***
556 :
557 : !!****f* m_raman/defwrite_raman_terms
558 : !! NAME
559 : !! defwrite_raman_terms
560 : !!
561 : !! FUNCTION
562 : !! Write the Raman susceptiblities for q=0 and also the phonon frequncies at gamma.
563 : !!
564 : !! INPUTS
565 : !! ncid=netcdf file id.
566 : !! rsus(3*natom,3,3)=List of Raman susceptibilities.
567 : !! natom=Number of atoms
568 : !!
569 : !! OUTPUT
570 : !! Only writing.
571 : !!
572 : !! SOURCE
573 :
574 2 : subroutine defwrite_raman_terms(ncid, natom, rsus, phfrq)
575 :
576 : !Arguments ------------------------------------
577 : !scalars
578 : integer,intent(in) :: ncid,natom
579 : !arrays
580 : real(dp),intent(in) :: rsus(3*natom,3,3)
581 : real(dp),intent(in) :: phfrq(3*natom)
582 :
583 : !Local variables-------------------------------
584 : !scalars
585 : integer :: ncerr, raman_sus_varid, phmodes_varid
586 : ! *************************************************************************
587 :
588 2 : NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.))
589 : ncerr = nctk_def_arrays(ncid, [ nctkarr_t("raman_sus", "dp", &
590 : "number_of_phonon_modes,number_of_cartesian_directions,number_of_cartesian_directions"), &
591 6 : nctkarr_t("gamma_phonon_modes", "dp", "number_of_phonon_modes")])
592 2 : NCF_CHECK(ncerr)
593 :
594 2 : NCF_CHECK(nctk_set_datamode(ncid))
595 :
596 2 : NCF_CHECK(nf90_inq_varid(ncid, "raman_sus", raman_sus_varid))
597 2 : NCF_CHECK(nf90_put_var(ncid,raman_sus_varid,rsus))
598 2 : NCF_CHECK(nf90_inq_varid(ncid, "gamma_phonon_modes", phmodes_varid))
599 14 : NCF_CHECK(nf90_put_var(ncid,phmodes_varid,phfrq*Ha_eV))
600 :
601 2 : end subroutine defwrite_raman_terms
602 : !!***
603 :
604 : end module m_raman
605 : !!***
|
