Branch data Line data Source code
1 : : !!****m* ABINIT/m_strain
2 : : !!
3 : : !! NAME
4 : : !! m_strain
5 : : !!
6 : : !! FUNCTION
7 : : !! Module for get the strain
8 : : !! Container type is defined
9 : : !!
10 : : !! COPYRIGHT
11 : : !! Copyright (C) 2010-2026 ABINIT group (AM)
12 : : !! This file is distributed under the terms of the
13 : : !! GNU General Public Licence, see ~abinit/COPYING
14 : : !! or http://www.gnu.org/copyleft/gpl.txt .
15 : : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
16 : : !!
17 : : !! SOURCE
18 : :
19 : : #if defined HAVE_CONFIG_H
20 : : #include "config.h"
21 : : #endif
22 : :
23 : : #include "abi_common.h"
24 : :
25 : : module m_strain
26 : :
27 : : use defs_basis
28 : : use m_errors
29 : : use m_abicore
30 : : use m_xmpi
31 : :
32 : : use m_matrix, only : matr3inv
33 : :
34 : : implicit none
35 : :
36 : : private :: strain_def2strain
37 : : private :: strain_strain2def
38 : : public :: strain_print
39 : : public :: strain_get
40 : : public :: strain_init
41 : : public :: strain_apply
42 : : !!***
43 : :
44 : : !!****t* m_strain/strain_type
45 : : !! NAME
46 : : !! strain_type
47 : : !!
48 : : !! FUNCTION
49 : : !! structure for a effective potential constructed.
50 : : !!
51 : : !! SOURCE
52 : :
53 : : type, public :: strain_type
54 : : character(len=fnlen) :: name
55 : : ! name of the strain (iso,uniaxial,shear...)
56 : :
57 : : real(dp) :: delta
58 : : ! Value of the strain
59 : :
60 : : integer :: direction
61 : : ! Direction of the strain (-1 if isostatic)
62 : :
63 : : real(dp) :: strain(3,3)
64 : : ! Matrix representing the strain
65 : :
66 : : end type strain_type
67 : : !!***
68 : :
69 : : CONTAINS !===========================================================================================
70 : :
71 : : !****f* m_strain/strain_init
72 : : !!
73 : : !! NAME
74 : : !! strain_init
75 : : !!
76 : : !! FUNCTION
77 : : !! routine to initialize strain structure
78 : : !!
79 : : !!
80 : : !! INPUTS
81 : : !! name = name of the perturbation
82 : : !! direction = direction of the perturbation
83 : : !! delta = delta to apply in the strain (in percent)
84 : : !!
85 : : !! OUTPUT
86 : : !! strain = structure with all information of strain
87 : : !!
88 : : !! SOURCE
89 : :
90 : 0 : subroutine strain_init(strain,delta,direction,name)
91 : :
92 : : !Arguments ------------------------------------
93 : : !scalars
94 : : character(len=fnlen),optional,intent(in) :: name
95 : : real(dp),optional,intent(in) :: delta
96 : : integer,optional,intent(in) :: direction
97 : : !array
98 : : type(strain_type),intent(out) :: strain
99 : : !Local variables-------------------------------
100 : : !scalar
101 : : !arrays
102 : : ! *************************************************************************
103 [ # # ]: 0 : if (present(name)) then
104 : 0 : strain%name = name
105 : : else
106 : 0 : strain%name = ''
107 : : end if
108 : :
109 [ # # ]: 0 : if (present(delta)) then
110 : 0 : strain%delta = delta
111 : : else
112 : 0 : strain%delta = zero
113 : : end if
114 : :
115 [ # # ]: 0 : if (present(direction)) then
116 : 0 : strain%direction = direction
117 : : else
118 : 0 : strain%direction = 0
119 : : end if
120 : :
121 : 0 : call strain_strain2def(strain%strain,strain)
122 : :
123 : 0 : end subroutine strain_init
124 : : !!***
125 : :
126 : : !****f* m_strain/strain_free
127 : : !!
128 : : !! NAME
129 : : !! strain_free
130 : : !!
131 : : !! FUNCTION
132 : : !! routine to free strain structure
133 : : !!
134 : : !!
135 : : !! INPUTS
136 : : !! strain = structure with all information of strain
137 : : !!
138 : : !! OUTPUT
139 : : !!
140 : : !! SOURCE
141 : :
142 : 0 : subroutine strain_free(strain)
143 : :
144 : : !Arguments ------------------------------------
145 : : !scalars
146 : : !array
147 : : type(strain_type),intent(inout) :: strain
148 : :
149 : : !Local variables-------------------------------
150 : : !scalar
151 : : !arrays
152 : : ! *************************************************************************
153 : :
154 : 0 : strain%name = ''
155 : 0 : strain%delta = zero
156 : 0 : strain%direction = 0
157 [ # # # # ]: 0 : strain%strain = zero
158 : :
159 : 0 : end subroutine strain_free
160 : : !!***
161 : :
162 : :
163 : : !****f* m_strain/strain_get
164 : : !!
165 : : !! NAME
166 : : !! strain_get
167 : : !!
168 : : !! FUNCTION
169 : : !! Get the strain for structure, compare to reference
170 : : !! structure and fill strain type
171 : : !!
172 : : !!
173 : : !! INPUTS
174 : : !! symmetrized = (optional) symmetrize the output
175 : : !!
176 : : !! OUTPUT
177 : : !! strain = structure with all information of strain
178 : : !!
179 : : !! SOURCE
180 : :
181 : 27999 : subroutine strain_get(strain,rprim,rprim_def,mat_delta,symmetrized)
182 : :
183 : : !Arguments ------------------------------------
184 : : !scalars
185 : : !array
186 : : type(strain_type),intent(inout) :: strain
187 : : real(dp),optional,intent(in) :: rprim(3,3),rprim_def(3,3), mat_delta(3,3)
188 : : logical,optional,intent(in) :: symmetrized
189 : : !Local variables-------------------------------
190 : : !scalar
191 : : integer :: i,j
192 : : logical :: symmetrized_in
193 : : character(len=500) :: message
194 : : !arrays
195 : : real(dp) :: mat_delta_tmp(3,3),rprim_inv(3,3)
196 : : real(dp) :: identity(3,3)
197 : : ! *************************************************************************
198 : :
199 : : !check inputs
200 : 27999 : symmetrized_in = .FALSE.
201 [ + + ]: 27999 : if(present(symmetrized)) then
202 : 13519 : symmetrized_in = symmetrized
203 : : end if
204 : :
205 : 27999 : if((present(rprim_def).and..not.present(rprim)).or.&
206 [ - + ]: 27999 : & (present(rprim).and..not.present(rprim_def))) then
207 : : write(message, '(a)' )&
208 : 0 : & ' strain_get: should give rprim_def and rprim as input of the routines'
209 : 0 : ABI_BUG(message)
210 : : end if
211 : :
212 [ + + ]: 27999 : if(present(rprim_def).and.present(rprim))then
213 : : mat_delta_tmp = zero
214 : : ! Fill the identity matrix
215 : 27038 : identity = zero
216 [ + + ]: 108152 : forall(i=1:3)identity(i,i)=1
217 : :
218 : 27038 : call matr3inv(rprim,rprim_inv)
219 [ + + + + : 1433014 : mat_delta_tmp = matmul(rprim_def,transpose(rprim_inv))-identity
+ + + + +
+ ]
220 : 27038 : identity = zero
221 [ + + ]: 108152 : do i=1,3
222 [ + + ]: 351494 : do j=1,3
223 [ + + ]: 324456 : if (abs(mat_delta_tmp(i,j))>tol10) then
224 : 76532 : identity(i,j) = mat_delta_tmp(i,j)
225 : : end if
226 : : end do
227 : : end do
228 : :
229 : 27038 : mat_delta_tmp = identity
230 : :
231 [ + - ]: 961 : else if (present(mat_delta)) then
232 : 961 : mat_delta_tmp = mat_delta
233 : :
234 : : else
235 : : write(message, '(a)' )&
236 : 0 : & ' strain_get: should give rprim_def or mat_delta as input of the routines'
237 : 0 : ABI_BUG(message)
238 : : end if
239 : :
240 [ - + ]: 27999 : if(symmetrized_in)then
241 : 0 : mat_delta_tmp(2,3) = (mat_delta_tmp(2,3) + mat_delta_tmp(3,2)) / 2
242 : 0 : mat_delta_tmp(3,1) = (mat_delta_tmp(3,1) + mat_delta_tmp(1,3)) / 2
243 : 0 : mat_delta_tmp(2,1) = (mat_delta_tmp(2,1) + mat_delta_tmp(1,2)) / 2
244 : :
245 : 0 : mat_delta_tmp(3,2) = mat_delta_tmp(2,3)
246 : 0 : mat_delta_tmp(1,3) = mat_delta_tmp(3,1)
247 : 0 : mat_delta_tmp(1,2) = mat_delta_tmp(2,1)
248 : :
249 : : end if
250 : :
251 : 27999 : call strain_def2strain(mat_delta_tmp,strain)
252 : :
253 : 27999 : end subroutine strain_get
254 : : !!***
255 : :
256 : : !****f* m_strain/strain_apply
257 : : !!
258 : : !! NAME
259 : : !! strain_get
260 : : !!
261 : : !! FUNCTION
262 : : !! Get the strain for structure, compare to reference
263 : : !! structure and fill strain type
264 : : !!
265 : : !!
266 : : !! INPUTS
267 : : !!
268 : : !! OUTPUT
269 : : !! strain = structure with all information of strain
270 : : !!
271 : : !! SOURCE
272 : :
273 : 0 : subroutine strain_apply(rprim,rprim_def,strain)
274 : :
275 : : !Arguments ------------------------------------
276 : : !scalars
277 : : !array
278 : : real(dp),intent(in) :: rprim(3,3)
279 : : real(dp),intent(out) :: rprim_def(3,3)
280 : : type(strain_type),intent(in) :: strain
281 : : !Local variables-------------------------------
282 : : !scalar
283 : : !integer :: i
284 : : !arrays
285 : : ! *************************************************************************
286 : :
287 : 0 : rprim_def(:,:) = zero
288 : : ! Fill the identity matrix
289 [ # # # # : 0 : rprim_def(:,:) = matmul(strain%strain(:,:),transpose(rprim(:,:)))
# # ]
290 : :
291 : 0 : end subroutine strain_apply
292 : : !!***
293 : :
294 : : !****f* m_strain/strain_def2strain
295 : : !!
296 : : !! NAME
297 : : !! strain_matdef2strain
298 : : !!
299 : : !! FUNCTION
300 : : !! transfer deformation matrix in structure strain
301 : : !!
302 : : !! INPUTS
303 : : !! rprim = contains
304 : : !!
305 : : !! OUTPUT
306 : : !!
307 : : !!
308 : : !! SOURCE
309 : :
310 : 27999 : subroutine strain_def2strain(mat_strain,strain)
311 : :
312 : : !Arguments ------------------------------------
313 : : !scalars
314 : : !array
315 : : real(dp),intent(in) :: mat_strain(3,3)
316 : : type(strain_type),intent(inout) :: strain
317 : : !Local variables-------------------------------
318 : : !scalar
319 : : !arrays
320 : : ! *************************************************************************
321 : 27999 : strain%name = ""
322 : 27999 : strain%delta = zero
323 : 27999 : strain%direction = 0
324 [ + + + + ]: 363987 : strain%strain = mat_strain
325 : :
326 [ + + + + : 241167 : if (all(abs(mat_strain)<tol10)) then
+ + + + ]
327 : 17764 : strain%name = "reference"
328 : : strain%delta = zero
329 : : strain%direction = 0
330 [ + + + + ]: 230932 : strain%strain = zero
331 : : else
332 : : if(abs(mat_strain(1,1))>tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))<tol10.and.&
333 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))>tol10.and.abs(mat_strain(2,3))<tol10.and.&
334 [ + - + + : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))>tol10) then
+ + + - +
- + - + -
+ - + - ]
335 [ + - + - ]: 1321 : if((mat_strain(1,1)-mat_strain(2,2))< tol10.and.&
336 : : & (mat_strain(1,1)-mat_strain(3,3))< tol10) then
337 : 1321 : strain%name = "isostatic"
338 : 1321 : strain%delta = mat_strain(1,1)
339 : 1321 : strain%direction = -1
340 [ + + + + ]: 17173 : strain%strain = mat_strain
341 : : end if
342 : : end if
343 : : if(abs(mat_strain(1,1))>tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))<tol10.and.&
344 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))<tol10.and.abs(mat_strain(2,3))<tol10.and.&
345 [ + - + + : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))<tol10) then
+ + + - -
+ - - - -
- - - - ]
346 : 0 : strain%name = "uniaxial"
347 : 0 : strain%delta = mat_strain(1,1)
348 : 0 : strain%direction = 1
349 [ # # # # ]: 0 : strain%strain = mat_strain
350 : : end if
351 : : if(abs(mat_strain(1,1))<tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))<tol10.and.&
352 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))>tol10.and.abs(mat_strain(2,3))<tol10.and.&
353 [ - + - - : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))<tol10) then
- - - - -
- - - - -
- - - - ]
354 : 0 : strain%name = "uniaxial"
355 : 0 : strain%delta = mat_strain(2,2)
356 : 0 : strain%direction = 2
357 [ # # # # ]: 0 : strain%strain = mat_strain
358 : : end if
359 : : if(abs(mat_strain(1,1))<tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))<tol10.and.&
360 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))<tol10.and.abs(mat_strain(2,3))<tol10.and.&
361 [ - + - - : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))>tol10) then
- - - - -
- - - - -
- - - - ]
362 : 0 : strain%name = "uniaxial"
363 : 0 : strain%delta = mat_strain(3,3)
364 : 0 : strain%direction = 3
365 [ # # # # ]: 0 : strain%strain = mat_strain
366 : : end if
367 : : if(abs(mat_strain(1,1))<tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))<tol10.and.&
368 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))<tol10.and.abs(mat_strain(2,3))>tol10.and.&
369 [ - + - - : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))>tol10.and.abs(mat_strain(3,3))<tol10) then
- - - - -
- - - - -
- - - - ]
370 [ # # ]: 0 : if (abs(mat_strain(3,2)-mat_strain(3,2))<tol10) then
371 : 0 : strain%name = "shear"
372 : 0 : strain%delta = mat_strain(3,2) * 2
373 : 0 : strain%direction = 4
374 [ # # # # ]: 0 : strain%strain = mat_strain
375 : : end if
376 : : end if
377 : : if(abs(mat_strain(1,1))<tol10.and.abs(mat_strain(1,2))<tol10.and.abs(mat_strain(1,3))>tol10.and.&
378 : : & abs(mat_strain(2,1))<tol10.and.abs(mat_strain(2,2))<tol10.and.abs(mat_strain(2,3))<tol10.and.&
379 [ - + - - : 10235 : & abs(mat_strain(3,1))>tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))<tol10) then
- - - - -
- - - - -
- - - - ]
380 [ # # ]: 0 : if (abs(mat_strain(3,1)-mat_strain(1,3))<tol10) then
381 : 0 : strain%name = "shear"
382 : 0 : strain%delta = mat_strain(3,1) * 2
383 : 0 : strain%direction = 5
384 [ # # # # ]: 0 : strain%strain = mat_strain
385 : : end if
386 : : end if
387 : : if(abs(mat_strain(1,1))<tol10.and.abs(mat_strain(1,2))>tol10.and.abs(mat_strain(1,3))<tol10.and.&
388 : : & abs(mat_strain(2,1))>tol10.and.abs(mat_strain(2,2))<tol10.and.abs(mat_strain(2,3))<tol10.and.&
389 [ - + - - : 10235 : & abs(mat_strain(3,1))<tol10.and.abs(mat_strain(3,2))<tol10.and.abs(mat_strain(3,3))<tol10) then
- - - - -
- - - - -
- - - - ]
390 [ # # ]: 0 : if (abs(mat_strain(1,2)-mat_strain(2,1))<tol10) then
391 : 0 : strain%name = "shear"
392 : 0 : strain%delta = mat_strain(2,1) * 2
393 : 0 : strain%direction = 6
394 [ # # # # ]: 0 : strain%strain = mat_strain
395 : : end if
396 : : end if
397 : : end if
398 : :
399 : 27999 : end subroutine strain_def2strain
400 : : !!***
401 : :
402 : : !****f* m_strain/strain_strain2def
403 : : !!
404 : : !! NAME
405 : : !! strain_matdef2strain
406 : : !!
407 : : !! FUNCTION
408 : : !! transfer deformation matrix in structure strain
409 : : !!
410 : : !! INPUTS
411 : : !! rprim = contains
412 : : !!
413 : : !! OUTPUT
414 : : !!
415 : : !!
416 : : !! SOURCE
417 : :
418 : 0 : subroutine strain_strain2def(mat_strain,strain)
419 : :
420 : : !Arguments ------------------------------------
421 : : !scalars
422 : : !array
423 : : real(dp),intent(out) :: mat_strain(3,3)
424 : : type(strain_type),intent(in) :: strain
425 : : !Local variables-------------------------------
426 : : !scalar
427 : : integer :: i
428 : : !arrays
429 : : ! *************************************************************************
430 : :
431 : 0 : mat_strain(:,:) = zero
432 [ # # ]: 0 : forall(i=1:3)mat_strain(i,i)=1
433 : :
434 : 0 : select case(strain%direction)
435 : : case(1)
436 : 0 : mat_strain(1,1) = mat_strain(1,1) + strain%delta
437 : : case(2)
438 : 0 : mat_strain(2,2) = mat_strain(2,2) + strain%delta
439 : : case(3)
440 : 0 : mat_strain(3,3) = mat_strain(3,3) + strain%delta
441 : : case(4)
442 : 0 : mat_strain(3,2) = strain%delta / 2
443 : 0 : mat_strain(2,3) = strain%delta / 2
444 : : case(5)
445 : 0 : mat_strain(3,1) = strain%delta / 2
446 : 0 : mat_strain(1,3) = strain%delta / 2
447 : : case(6)
448 : 0 : mat_strain(2,1) = strain%delta / 2
449 [ # # # # : 0 : mat_strain(1,2) = strain%delta / 2
# # # ]
450 : : end select
451 : :
452 : 0 : end subroutine strain_strain2def
453 : : !!***
454 : :
455 : : !****f* m_strain/strain_print
456 : : !!
457 : : !! NAME
458 : : !! strain_print
459 : : !!
460 : : !! FUNCTION
461 : : !! print the structure strain
462 : : !!
463 : : !! INPUTS
464 : : !!
465 : : !! OUTPUT
466 : : !! eff_pot = supercell structure with data to be output
467 : : !!
468 : : !! SOURCE
469 : :
470 : 2056 : subroutine strain_print(strain)
471 : :
472 : : !Arguments ------------------------------------
473 : : !scalars
474 : : !array
475 : : type(strain_type),intent(in) :: strain
476 : : !Local variables-------------------------------
477 : : !scalar
478 : : integer :: ii
479 : : character(len=500) :: message
480 : : !arrays
481 : : ! *************************************************************************
482 : :
483 [ + + ]: 2056 : if(strain%name == "reference") then
484 : 1081 : write(message,'(4a)') ch10,' no strain found:',&
485 : 2162 : & ' This structure is equivalent to the reference structure',ch10
486 : 1081 : call wrtout(std_out,message,'COLL')
487 : : else
488 [ + + ]: 975 : if(strain%name /= "") then
489 : : write(message,'(3a,I2,a,(ES10.2),a)') &
490 : 20 : & ' The strain is ',trim(strain%name),' type in the direction ',&
491 : 40 : & strain%direction,' with delta of ',strain%delta, ':'
492 : 20 : call wrtout(std_out,message,'COLL')
493 : 20 : call wrtout(ab_out,message,'COLL')
494 [ + + ]: 80 : do ii = 1,3
495 : 60 : write(message,'(3es17.8)') strain%strain(ii,1),strain%strain(ii,2),strain%strain(ii,3)
496 : 80 : call wrtout(std_out,message,'COLL')
497 : : end do
498 : : else
499 : 955 : write(message,'(a)') ' Strain does not correspond to standard strain:'
500 : 955 : call wrtout(std_out,message,'COLL')
501 [ + + ]: 3820 : do ii = 1,3
502 : 2865 : write(message,'(3es17.8)') strain%strain(ii,1),strain%strain(ii,2),strain%strain(ii,3)
503 : 3820 : call wrtout(std_out,message,'COLL')
504 : : end do
505 : : end if
506 : : end if
507 : 2056 : end subroutine strain_print
508 : : !!***
509 : :
510 : 0 : end module m_strain
511 : :
|
