Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief unit conversion facility
10 : !>
11 : !> Units are complex, this module does not try to be very smart, for
12 : !> example SI prefixes are not supported automatically, and
13 : !> which kinds are really basic can change depending on the system of
14 : !> units chosen, and equivalences are not always catched.
15 : !>
16 : !> This is thought as a simple conversion facility for the input and output.
17 : !> If you need something more you are probably better off using the
18 : !> physcon module directly.
19 : !> \note
20 : !> One design choice was not to use dynamically allocated elements to
21 : !> reduce the possibility of leaks.
22 : !> Needs to be extended (for example charge, dipole,...)
23 : !> I just added the units and kinds that I needed.
24 : !> Used by the parser
25 : !> Should keep an unsorted/uncompressed version for nicer labels?
26 : !> \par History
27 : !> 01.2005 created [fawzi]
28 : !> \author fawzi
29 : ! **************************************************************************************************
30 : MODULE cp_units
31 :
32 : USE cp_log_handling, ONLY: cp_to_string
33 : USE kinds, ONLY: dp
34 : USE mathconstants, ONLY: radians,&
35 : twopi
36 : USE physcon, ONLY: &
37 : atm, bar, bohr, e_mass, evolt, femtoseconds, joule, kcalmol, kelvin, kjmol, massunit, &
38 : newton, pascal, picoseconds, seconds, wavenumbers
39 : USE string_utilities, ONLY: compress,&
40 : s2a,&
41 : uppercase
42 : #include "../base/base_uses.f90"
43 :
44 : IMPLICIT NONE
45 : PRIVATE
46 :
47 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
48 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_units'
49 :
50 : INTEGER, PARAMETER, PUBLIC :: cp_ukind_none = 0, &
51 : cp_ukind_energy = 1, &
52 : cp_ukind_length = 2, &
53 : cp_ukind_temperature = 3, &
54 : cp_ukind_angle = 4, &
55 : cp_ukind_pressure = 5, &
56 : cp_ukind_time = 6, &
57 : cp_ukind_mass = 7, &
58 : cp_ukind_undef = 8, &
59 : cp_ukind_potential = 9, &
60 : cp_ukind_force = 10, &
61 : cp_ukind_efield = 11, &
62 : cp_ukind_max = 11
63 :
64 : ! General
65 : INTEGER, PARAMETER, PUBLIC :: cp_units_none = 100, &
66 : cp_units_au = 101
67 : ! Mass
68 : INTEGER, PARAMETER, PUBLIC :: cp_units_m_e = 110, &
69 : cp_units_amu = 111, &
70 : cp_units_kg = 112
71 : ! Energy
72 : INTEGER, PARAMETER, PUBLIC :: cp_units_hartree = 130, &
73 : cp_units_wavenum = 131, &
74 : cp_units_joule = 132, &
75 : cp_units_kcalmol = 133, &
76 : cp_units_Ry = 134, &
77 : cp_units_eV = 135, &
78 : cp_units_kjmol = 136, &
79 : cp_units_jmol = 137, &
80 : cp_units_keV = 138
81 :
82 : ! Length
83 : INTEGER, PARAMETER, PUBLIC :: cp_units_bohr = 140, &
84 : cp_units_angstrom = 141, &
85 : cp_units_m = 142, &
86 : cp_units_pm = 143, &
87 : cp_units_nm = 144
88 :
89 : ! Temperature
90 : INTEGER, PARAMETER, PUBLIC :: cp_units_k = 150
91 :
92 : ! Pressure
93 : INTEGER, PARAMETER, PUBLIC :: cp_units_bar = 161
94 : INTEGER, PARAMETER, PUBLIC :: cp_units_atm = 162
95 : INTEGER, PARAMETER, PUBLIC :: cp_units_kbar = 163
96 : INTEGER, PARAMETER, PUBLIC :: cp_units_Pa = 164
97 : INTEGER, PARAMETER, PUBLIC :: cp_units_MPa = 165
98 : INTEGER, PARAMETER, PUBLIC :: cp_units_GPa = 166
99 :
100 : ! Angles
101 : INTEGER, PARAMETER, PUBLIC :: cp_units_rad = 170, &
102 : cp_units_deg = 171
103 :
104 : ! Time
105 : INTEGER, PARAMETER, PUBLIC :: cp_units_fs = 180, &
106 : cp_units_s = 181, &
107 : cp_units_wn = 182, &
108 : cp_units_ps = 183
109 :
110 : ! Potential
111 : INTEGER, PARAMETER, PUBLIC :: cp_units_volt = 190
112 :
113 : ! Force
114 : INTEGER, PARAMETER, PUBLIC :: cp_units_Newton = 200, &
115 : cp_units_mNewton = 201
116 :
117 : ! Electric Field
118 : INTEGER, PARAMETER, PUBLIC :: cp_units_volt_per_m = 202, &
119 : cp_units_volt_per_nm = 203, &
120 : cp_units_volt_per_angstrom = 204
121 :
122 : INTEGER, PARAMETER, PUBLIC :: cp_unit_max_kinds = 8, cp_unit_basic_desc_length = 15, &
123 : cp_unit_desc_length = cp_unit_max_kinds*cp_unit_basic_desc_length
124 :
125 : PUBLIC :: cp_unit_type, cp_unit_set_type
126 : PUBLIC :: cp_unit_create, cp_unit_release, &
127 : cp_unit_to_cp2k, cp_unit_from_cp2k, cp_unit_desc, &
128 : cp_unit_set_create, cp_unit_set_release, &
129 : cp_unit_to_cp2k1, cp_unit_from_cp2k1, cp_unit_compatible, export_units_as_xml
130 :
131 : ! **************************************************************************************************
132 : !> \brief stores a unit
133 : !> \param kind the kind of unit (energy, length,...)
134 : !> \param unit the actual unit (Joule, eV,...)
135 : !> \author fawzi
136 : ! **************************************************************************************************
137 : TYPE cp_unit_type
138 : INTEGER :: n_kinds = -1
139 : INTEGER, DIMENSION(cp_unit_max_kinds):: kind_id = -1, unit_id = -1, power = -1
140 : END TYPE cp_unit_type
141 :
142 : ! **************************************************************************************************
143 : !> \brief represent a pointer to a unit (to build arrays of pointers)
144 : !> \param unit the pointer to the unit
145 : !> \author fawzi
146 : ! **************************************************************************************************
147 : TYPE cp_unit_p_type
148 : TYPE(cp_unit_type), POINTER :: unit => NULL()
149 : END TYPE cp_unit_p_type
150 :
151 : ! **************************************************************************************************
152 : !> \brief stores the default units to be used
153 : !> \author fawzi
154 : ! **************************************************************************************************
155 : TYPE cp_unit_set_type
156 : TYPE(cp_unit_p_type), DIMENSION(cp_ukind_max) :: units = cp_unit_p_type()
157 : END TYPE cp_unit_set_type
158 :
159 : CONTAINS
160 :
161 : ! **************************************************************************************************
162 : !> \brief creates a unit parsing a string
163 : !> \param unit the unit to initialize
164 : !> \param string the string containing the description of the unit
165 : !> \author fawzi
166 : ! **************************************************************************************************
167 193289850 : SUBROUTINE cp_unit_create(unit, string)
168 : TYPE(cp_unit_type), INTENT(OUT) :: unit
169 : CHARACTER(len=*), INTENT(in) :: string
170 :
171 : CHARACTER(LEN=40) :: formatstr
172 : CHARACTER(LEN=cp_unit_desc_length) :: desc
173 7731594 : CHARACTER(LEN=LEN(string)) :: unit_string
174 : INTEGER :: i_high, i_low, i_unit, len_string, &
175 : next_power
176 : INTEGER, DIMENSION(cp_unit_max_kinds) :: kind_id, power, unit_id
177 :
178 69584346 : unit_id = cp_units_none
179 7731594 : kind_id = cp_ukind_none
180 7731594 : power = 0
181 7731594 : i_low = 1
182 7731594 : i_high = 1
183 7731594 : len_string = LEN(string)
184 7731594 : i_unit = 0
185 7731594 : next_power = 1
186 7731594 : DO WHILE (i_low < len_string)
187 7274440 : IF (string(i_low:i_low) /= ' ') EXIT
188 7731594 : i_low = i_low + 1
189 : END DO
190 : i_high = i_low
191 51275362 : DO WHILE (i_high <= len_string)
192 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
193 44823462 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
194 49995723 : i_high = i_high + 1
195 : END DO
196 : DO
197 8362192 : IF (i_high <= i_low .OR. i_low > len_string) EXIT
198 8354027 : i_unit = i_unit + 1
199 8354027 : IF (i_unit > cp_unit_max_kinds) THEN
200 0 : CPABORT("Maximum number of combined units exceeded")
201 0 : EXIT
202 : END IF
203 : ! read unit
204 8354027 : unit_string = string(i_low:i_high - 1)
205 8354027 : CALL uppercase(unit_string)
206 8665351 : SELECT CASE (TRIM(unit_string))
207 : CASE ("INTERNAL_CP2K")
208 311324 : unit_id(i_unit) = cp_units_none
209 311324 : kind_id(i_unit) = cp_ukind_undef
210 : CASE ("HARTREE")
211 293826 : unit_id(i_unit) = cp_units_hartree
212 293826 : kind_id(i_unit) = cp_ukind_energy
213 : CASE ("AU_E")
214 36790 : unit_id(i_unit) = cp_units_au
215 36790 : kind_id(i_unit) = cp_ukind_energy
216 : CASE ("WAVENUMBER_E")
217 0 : unit_id(i_unit) = cp_units_wavenum
218 0 : kind_id(i_unit) = cp_ukind_energy
219 : CASE ("JOULE", "J")
220 0 : unit_id(i_unit) = cp_units_joule
221 0 : kind_id(i_unit) = cp_ukind_energy
222 : CASE ("KCALMOL")
223 234924 : unit_id(i_unit) = cp_units_kcalmol
224 234924 : kind_id(i_unit) = cp_ukind_energy
225 : CASE ("KJMOL")
226 2248 : unit_id(i_unit) = cp_units_kjmol
227 2248 : kind_id(i_unit) = cp_ukind_energy
228 : CASE ("JMOL")
229 0 : unit_id(i_unit) = cp_units_jmol
230 0 : kind_id(i_unit) = cp_ukind_energy
231 : CASE ("RY")
232 51565 : unit_id(i_unit) = cp_units_Ry
233 51565 : kind_id(i_unit) = cp_ukind_energy
234 : CASE ("EV")
235 681424 : unit_id(i_unit) = cp_units_eV
236 681424 : kind_id(i_unit) = cp_ukind_energy
237 : CASE ("KEV")
238 21700 : unit_id(i_unit) = cp_units_keV
239 21700 : kind_id(i_unit) = cp_ukind_energy
240 : CASE ("K_E")
241 98696 : unit_id(i_unit) = cp_units_k
242 98696 : kind_id(i_unit) = cp_ukind_energy
243 : CASE ("ENERGY")
244 0 : unit_id(i_unit) = cp_units_none
245 0 : kind_id(i_unit) = cp_ukind_energy
246 : CASE ("AU_L")
247 278 : unit_id(i_unit) = cp_units_au
248 278 : kind_id(i_unit) = cp_ukind_length
249 : CASE ("BOHR")
250 607855 : unit_id(i_unit) = cp_units_bohr
251 607855 : kind_id(i_unit) = cp_ukind_length
252 : CASE ("M")
253 24497 : unit_id(i_unit) = cp_units_m
254 24497 : kind_id(i_unit) = cp_ukind_length
255 : CASE ("PM")
256 2 : unit_id(i_unit) = cp_units_pm
257 2 : kind_id(i_unit) = cp_ukind_length
258 : CASE ("NM")
259 11417 : unit_id(i_unit) = cp_units_nm
260 11417 : kind_id(i_unit) = cp_ukind_length
261 : CASE ("ANGSTROM")
262 3602252 : unit_id(i_unit) = cp_units_angstrom
263 3602252 : kind_id(i_unit) = cp_ukind_length
264 : CASE ("LENGTH")
265 0 : unit_id(i_unit) = cp_units_none
266 0 : kind_id(i_unit) = cp_ukind_length
267 : CASE ("K", "K_TEMP")
268 484822 : unit_id(i_unit) = cp_units_k
269 484822 : kind_id(i_unit) = cp_ukind_temperature
270 : CASE ("AU_TEMP")
271 0 : unit_id(i_unit) = cp_units_au
272 0 : kind_id(i_unit) = cp_ukind_temperature
273 : CASE ("TEMPERATURE")
274 0 : unit_id(i_unit) = cp_units_none
275 0 : kind_id(i_unit) = cp_ukind_temperature
276 : CASE ("ATM")
277 0 : unit_id(i_unit) = cp_units_atm
278 0 : kind_id(i_unit) = cp_ukind_pressure
279 : CASE ("BAR")
280 83711 : unit_id(i_unit) = cp_units_bar
281 83711 : kind_id(i_unit) = cp_ukind_pressure
282 : CASE ("KBAR")
283 18 : unit_id(i_unit) = cp_units_kbar
284 18 : kind_id(i_unit) = cp_ukind_pressure
285 : CASE ("PA")
286 2716 : unit_id(i_unit) = cp_units_Pa
287 2716 : kind_id(i_unit) = cp_ukind_pressure
288 : CASE ("MPA")
289 0 : unit_id(i_unit) = cp_units_MPa
290 0 : kind_id(i_unit) = cp_ukind_pressure
291 : CASE ("GPA")
292 1666 : unit_id(i_unit) = cp_units_GPa
293 1666 : kind_id(i_unit) = cp_ukind_pressure
294 : CASE ("AU_P")
295 0 : unit_id(i_unit) = cp_units_au
296 0 : kind_id(i_unit) = cp_ukind_pressure
297 : CASE ("PRESSURE")
298 0 : unit_id(i_unit) = cp_units_none
299 0 : kind_id(i_unit) = cp_ukind_pressure
300 : CASE ("RAD")
301 87249 : unit_id(i_unit) = cp_units_rad
302 87249 : kind_id(i_unit) = cp_ukind_angle
303 : CASE ("DEG")
304 151576 : unit_id(i_unit) = cp_units_deg
305 151576 : kind_id(i_unit) = cp_ukind_angle
306 : CASE ("ANGLE")
307 0 : unit_id(i_unit) = cp_units_none
308 0 : kind_id(i_unit) = cp_ukind_angle
309 : CASE ("S")
310 38071 : unit_id(i_unit) = cp_units_s
311 38071 : kind_id(i_unit) = cp_ukind_time
312 : CASE ("FS")
313 1021278 : unit_id(i_unit) = cp_units_fs
314 1021278 : kind_id(i_unit) = cp_ukind_time
315 : CASE ("PS")
316 365474 : unit_id(i_unit) = cp_units_ps
317 365474 : kind_id(i_unit) = cp_ukind_time
318 : CASE ("WAVENUMBER_T")
319 34 : unit_id(i_unit) = cp_units_wn
320 34 : kind_id(i_unit) = cp_ukind_time
321 : CASE ("AU_T")
322 92131 : unit_id(i_unit) = cp_units_au
323 92131 : kind_id(i_unit) = cp_ukind_time
324 : CASE ("TIME")
325 0 : unit_id(i_unit) = cp_units_none
326 0 : kind_id(i_unit) = cp_ukind_time
327 : CASE ("KG")
328 0 : unit_id(i_unit) = cp_units_kg
329 0 : kind_id(i_unit) = cp_ukind_mass
330 : CASE ("AMU")
331 1514 : unit_id(i_unit) = cp_units_amu
332 1514 : kind_id(i_unit) = cp_ukind_mass
333 : CASE ("M_E")
334 0 : unit_id(i_unit) = cp_units_m_e
335 0 : kind_id(i_unit) = cp_ukind_mass
336 : CASE ("AU_M")
337 21696 : unit_id(i_unit) = cp_units_au
338 21696 : kind_id(i_unit) = cp_ukind_mass
339 : CASE ("MASS")
340 0 : unit_id(i_unit) = cp_units_none
341 0 : kind_id(i_unit) = cp_ukind_mass
342 : CASE ("VOLT")
343 17750 : unit_id(i_unit) = cp_units_volt
344 17750 : kind_id(i_unit) = cp_ukind_potential
345 : CASE ("AU_POT")
346 0 : unit_id(i_unit) = cp_units_au
347 0 : kind_id(i_unit) = cp_ukind_potential
348 : CASE ("POTENTIAL")
349 0 : unit_id(i_unit) = cp_units_none
350 0 : kind_id(i_unit) = cp_ukind_potential
351 : CASE ("N", "NEWTON")
352 12 : unit_id(i_unit) = cp_units_Newton
353 12 : kind_id(i_unit) = cp_ukind_force
354 : CASE ("MN", "MNEWTON")
355 2776 : unit_id(i_unit) = cp_units_mNewton
356 2776 : kind_id(i_unit) = cp_ukind_force
357 : CASE ("AU_F")
358 0 : unit_id(i_unit) = cp_units_au
359 0 : kind_id(i_unit) = cp_ukind_force
360 : CASE ("FORCE")
361 0 : unit_id(i_unit) = cp_units_none
362 0 : kind_id(i_unit) = cp_ukind_force
363 : CASE ("VM-1", "VOLT_PER_M")
364 2735 : unit_id(i_unit) = cp_units_volt_per_m
365 2735 : kind_id(i_unit) = cp_ukind_efield
366 : CASE ("VNM-1", "VOLT_PER_NM")
367 0 : unit_id(i_unit) = cp_units_volt_per_nm
368 0 : kind_id(i_unit) = cp_ukind_efield
369 : CASE ("VA-1", "VOLT_PER_ANGSTROM")
370 0 : unit_id(i_unit) = cp_units_volt_per_angstrom
371 0 : kind_id(i_unit) = cp_ukind_efield
372 : CASE ("AU_EFIELD")
373 0 : unit_id(i_unit) = cp_units_au
374 0 : kind_id(i_unit) = cp_ukind_efield
375 : CASE ("EFIELD")
376 0 : unit_id(i_unit) = cp_units_none
377 0 : kind_id(i_unit) = cp_ukind_efield
378 : CASE ("AU")
379 : CALL cp_abort(__LOCATION__, &
380 : "au unit without specifying its kind not accepted, use "// &
381 : "(au_e, au_f, au_t, au_temp, au_l, au_m, au_p, au_pot, "// &
382 0 : "au_efield)")
383 : CASE default
384 8354027 : CPABORT("Unknown unit: "//string(i_low:i_high - 1))
385 : END SELECT
386 8354027 : power(i_unit) = next_power
387 : ! parse op
388 8354027 : i_low = i_high
389 8377041 : DO WHILE (i_low <= len_string)
390 1893000 : IF (string(i_low:i_low) /= ' ') EXIT
391 8377041 : i_low = i_low + 1
392 : END DO
393 : i_high = i_low
394 8354027 : DO WHILE (i_high <= len_string)
395 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
396 1869986 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
397 0 : i_high = i_high + 1
398 : END DO
399 8354027 : IF (i_high < i_low .OR. i_low > len_string) EXIT
400 :
401 1869986 : IF (i_high <= len_string) THEN
402 1869986 : IF (string(i_low:i_high) == '^') THEN
403 1245057 : i_low = i_high + 1
404 1245057 : DO WHILE (i_low <= len_string)
405 1245057 : IF (string(i_low:i_low) /= ' ') EXIT
406 1245057 : i_low = i_low + 1
407 : END DO
408 : i_high = i_low
409 3682786 : DO WHILE (i_high <= len_string)
410 1239388 : SELECT CASE (string(i_high:i_high))
411 : CASE ('+', '-', '0', '1', '2', '3', '4', '5', '6', '7', '8', '9')
412 2437729 : i_high = i_high + 1
413 : CASE default
414 2443398 : EXIT
415 : END SELECT
416 : END DO
417 1245057 : IF (i_high <= i_low .OR. i_low > len_string) THEN
418 0 : CPABORT("an integer number is expected after a '^'")
419 0 : EXIT
420 : END IF
421 1245057 : formatstr = "(i"//cp_to_string(i_high - i_low + 1)//")"
422 : READ (string(i_low:i_high - 1), formatstr) &
423 1245057 : next_power
424 1245057 : power(i_unit) = power(i_unit)*next_power
425 : ! next op
426 1245057 : i_low = i_high
427 1249977 : DO WHILE (i_low < len_string)
428 10572 : IF (string(i_low:i_low) /= ' ') EXIT
429 1249977 : i_low = i_low + 1
430 : END DO
431 : i_high = i_low
432 1246251 : DO WHILE (i_high <= len_string)
433 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
434 6661 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
435 1211 : i_high = i_high + 1
436 : END DO
437 : END IF
438 : END IF
439 1869986 : IF (i_low > len_string) EXIT
440 630598 : next_power = 1
441 8362192 : IF (i_high <= len_string) THEN
442 630396 : IF (string(i_low:i_high) == "*" .OR. string(i_low:i_high) == '/') THEN
443 630379 : IF (string(i_low:i_high) == '/') next_power = -1
444 630379 : i_low = i_high + 1
445 630379 : DO WHILE (i_low <= len_string)
446 630379 : IF (string(i_low:i_low) /= ' ') EXIT
447 630379 : i_low = i_low + 1
448 : END DO
449 : i_high = i_low
450 3118677 : DO WHILE (i_high <= len_string)
451 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
452 3078647 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
453 2528347 : i_high = i_high + 1
454 : END DO
455 : END IF
456 : END IF
457 : END DO
458 : CALL cp_unit_create2(unit, kind_id=kind_id, unit_id=unit_id, &
459 7731594 : power=power)
460 7731594 : desc = cp_unit_desc(unit)
461 7731594 : END SUBROUTINE cp_unit_create
462 :
463 : ! **************************************************************************************************
464 : !> \brief creates and initializes the given unit of mesure (performs some error
465 : !> check)
466 : !> \param unit the unit descriptor to be initialized
467 : !> \param kind_id the kind of unit (length,energy,...), use the constants
468 : !> cp_ukind_*
469 : !> \param unit_id the actual unit (use constants cp_units_*)
470 : !> \param power ...
471 : !> \author fawzi
472 : ! **************************************************************************************************
473 220121804 : SUBROUTINE cp_unit_create2(unit, kind_id, unit_id, power)
474 : TYPE(cp_unit_type), INTENT(OUT) :: unit
475 : INTEGER, DIMENSION(:), INTENT(in) :: kind_id, unit_id
476 : INTEGER, DIMENSION(:), INTENT(in), OPTIONAL :: power
477 :
478 : INTEGER :: i, j, max_kind, max_pos
479 : LOGICAL :: repeat
480 :
481 7861493 : CPASSERT(SIZE(kind_id) <= cp_unit_max_kinds)
482 7861493 : CPASSERT(SIZE(unit_id) <= cp_unit_max_kinds)
483 69844144 : unit%kind_id(1:SIZE(kind_id)) = kind_id
484 8770786 : unit%kind_id(SIZE(kind_id) + 1:) = cp_ukind_none
485 69844144 : unit%unit_id(1:SIZE(unit_id)) = unit_id
486 16632279 : unit%unit_id(SIZE(unit_id):) = cp_units_none
487 7861493 : IF (PRESENT(power)) THEN
488 69844144 : unit%power(1:SIZE(power)) = power
489 8770786 : unit%power(SIZE(power) + 1:) = 0
490 70753437 : DO i = 1, SIZE(unit%power)
491 70753437 : IF (unit%power(i) == 0) THEN
492 54408018 : unit%kind_id(i) = cp_ukind_none
493 54408018 : unit%unit_id(i) = cp_units_none
494 : END IF
495 : END DO
496 : ELSE
497 0 : DO i = 1, SIZE(unit%power)
498 0 : IF (unit%unit_id(i) /= 0) THEN
499 0 : unit%power(i) = 1
500 : ELSE
501 0 : unit%power(i) = 0
502 : END IF
503 : END DO
504 : END IF
505 :
506 : ! remove unnecessary units
507 : ! reorder & compress
508 7861493 : unit%n_kinds = 0
509 70753437 : DO i = 1, SIZE(unit%kind_id)
510 : ! find max and compress in the rest
511 : DO
512 62891944 : max_kind = unit%kind_id(i)
513 62891944 : max_pos = i
514 62891944 : repeat = .FALSE.
515 283013748 : DO j = i + 1, SIZE(unit%kind_id)
516 283013748 : IF (unit%kind_id(j) >= max_kind) THEN
517 161992474 : IF (unit%kind_id(j) /= 0 .AND. unit%kind_id(j) == max_kind .AND. &
518 : unit%unit_id(j) == unit%unit_id(max_pos)) THEN
519 0 : unit%power(max_pos) = unit%power(max_pos) + unit%power(j)
520 0 : unit%kind_id(j) = cp_ukind_none
521 0 : unit%unit_id(j) = cp_units_none
522 0 : unit%power(j) = 0
523 0 : IF (unit%power(max_pos) == 0) THEN
524 0 : unit%kind_id(max_pos) = cp_ukind_none
525 0 : unit%unit_id(max_pos) = cp_units_none
526 0 : unit%power(max_pos) = 0
527 0 : repeat = .TRUE.
528 0 : EXIT
529 : END IF
530 161992474 : ELSE IF (unit%kind_id(j) > max_kind .OR. &
531 : (unit%kind_id(j) == max_kind .AND. &
532 : unit%unit_id(j) > unit%unit_id(max_pos))) THEN
533 624851 : max_kind = unit%kind_id(j)
534 624851 : max_pos = j
535 : END IF
536 : END IF
537 : END DO
538 62891944 : IF (.NOT. repeat) EXIT
539 : END DO
540 62891944 : IF (max_kind /= 0) unit%n_kinds = unit%n_kinds + 1
541 : ! put the max at pos i
542 62891944 : IF (max_pos /= i) THEN
543 622131 : unit%kind_id(max_pos) = unit%kind_id(i)
544 622131 : unit%kind_id(i) = max_kind
545 622131 : max_kind = unit%unit_id(max_pos)
546 622131 : unit%unit_id(max_pos) = unit%unit_id(i)
547 622131 : unit%unit_id(i) = max_kind
548 622131 : max_kind = unit%power(max_pos)
549 622131 : unit%power(max_pos) = unit%power(i)
550 622131 : unit%power(i) = max_kind
551 : END IF
552 : ! check unit
553 : CALL cp_basic_unit_check(basic_kind=unit%kind_id(i), &
554 70753437 : basic_unit=unit%unit_id(i))
555 : END DO
556 7861493 : END SUBROUTINE cp_unit_create2
557 :
558 : ! **************************************************************************************************
559 : !> \brief releases the given unit
560 : !> \param unit the unit to release
561 : !> \author fawzi
562 : !> \note
563 : !> at the moment not needed, there for completeness
564 : ! **************************************************************************************************
565 7861493 : ELEMENTAL SUBROUTINE cp_unit_release(unit)
566 : TYPE(cp_unit_type), INTENT(IN) :: unit
567 :
568 : MARK_USED(unit)
569 :
570 7861493 : END SUBROUTINE cp_unit_release
571 :
572 : ! **************************************************************************************************
573 : !> \brief controls that the kind and contains meaningful information
574 : !> \param basic_kind the kind of the unit
575 : !> \param basic_unit the unit to check
576 : !> \author fawzi
577 : ! **************************************************************************************************
578 62891944 : SUBROUTINE cp_basic_unit_check(basic_kind, basic_unit)
579 : INTEGER, INTENT(in) :: basic_kind, basic_unit
580 :
581 63215077 : SELECT CASE (basic_kind)
582 : CASE (cp_ukind_undef)
583 1756115 : SELECT CASE (basic_unit)
584 : CASE (cp_units_none)
585 : CASE default
586 323133 : CPABORT("unknown undef unit:"//TRIM(cp_to_string(basic_unit)))
587 : END SELECT
588 : CASE (cp_ukind_energy)
589 5691092 : SELECT CASE (basic_unit)
590 : CASE (cp_units_hartree, cp_units_wavenum, cp_units_joule, cp_units_kcalmol, &
591 : cp_units_kjmol, cp_units_Ry, cp_units_eV, cp_units_keV, cp_units_au, cp_units_k, &
592 : cp_units_jmol, cp_units_none)
593 : CASE default
594 1432982 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
595 : END SELECT
596 : CASE (cp_ukind_length)
597 4754741 : SELECT CASE (basic_unit)
598 : CASE (cp_units_bohr, cp_units_angstrom, cp_units_au, cp_units_none, cp_units_m, &
599 : cp_units_pm, cp_units_nm)
600 : CASE default
601 4258110 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
602 : END SELECT
603 : CASE (cp_ukind_temperature)
604 596551 : SELECT CASE (basic_unit)
605 : CASE (cp_units_k, cp_units_au, cp_units_none)
606 : CASE default
607 496631 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
608 : END SELECT
609 : CASE (cp_ukind_pressure)
610 350554 : SELECT CASE (basic_unit)
611 : CASE (cp_units_bar, cp_units_atm, cp_units_kbar, cp_units_Pa, cp_units_MPa, cp_units_GPa, cp_units_au, cp_units_none)
612 : CASE default
613 99920 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
614 : END SELECT
615 : CASE (cp_ukind_angle)
616 1779431 : SELECT CASE (basic_unit)
617 : CASE (cp_units_rad, cp_units_deg, cp_units_none)
618 : CASE default
619 250634 : CPABORT("unknown angle unit:"//TRIM(cp_to_string(basic_unit)))
620 : END SELECT
621 : CASE (cp_ukind_time)
622 1563816 : SELECT CASE (basic_unit)
623 : CASE (cp_units_s, cp_units_fs, cp_units_ps, cp_units_au, cp_units_wn, cp_units_none)
624 : CASE default
625 1528797 : CPABORT("unknown time unit:"//TRIM(cp_to_string(basic_unit)))
626 : END SELECT
627 : CASE (cp_ukind_mass)
628 64578 : SELECT CASE (basic_unit)
629 : CASE (cp_units_kg, cp_units_amu, cp_units_m_e, cp_units_au, cp_units_none)
630 : CASE default
631 35019 : CPABORT("unknown mass unit:"//TRIM(cp_to_string(basic_unit)))
632 : END SELECT
633 : CASE (cp_ukind_potential)
634 44156 : SELECT CASE (basic_unit)
635 : CASE (cp_units_volt, cp_units_au, cp_units_none)
636 : CASE default
637 29559 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
638 : END SELECT
639 : CASE (cp_ukind_force)
640 29141 : SELECT CASE (basic_unit)
641 : CASE (cp_units_Newton, cp_units_mNewton, cp_units_au, cp_units_none)
642 : CASE default
643 14597 : CPABORT("unknown force unit:"//TRIM(cp_to_string(basic_unit)))
644 : END SELECT
645 : CASE (cp_ukind_efield)
646 54422562 : SELECT CASE (basic_unit)
647 : CASE (cp_units_volt_per_m, cp_units_volt_per_nm, &
648 : cp_units_volt_per_angstrom, cp_units_au, cp_units_none)
649 : CASE default
650 14544 : CPABORT("unknown electric field unit:"//TRIM(cp_to_string(basic_unit)))
651 : END SELECT
652 : CASE (cp_ukind_none)
653 54408018 : IF (basic_unit /= cp_units_none) THEN
654 : CALL cp_abort(__LOCATION__, &
655 : "if the kind of the unit is none also unit must be undefined,not:" &
656 0 : //TRIM(cp_to_string(basic_unit)))
657 : END IF
658 : CASE default
659 62891944 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
660 : END SELECT
661 62891944 : END SUBROUTINE cp_basic_unit_check
662 :
663 : ! **************************************************************************************************
664 : !> \brief converts a value to the internal cp2k units
665 : !> \param value the value to convert
666 : !> \param basic_kind the kind of the unit of the value
667 : !> \param basic_unit the unit of the value
668 : !> \param power the power of the unit (defaults to 1)
669 : !> \return ...
670 : !> \author fawzi
671 : ! **************************************************************************************************
672 5208477 : FUNCTION cp_basic_unit_to_cp2k(value, basic_kind, basic_unit, power) RESULT(res)
673 : REAL(kind=dp), INTENT(in) :: value
674 : INTEGER, INTENT(in) :: basic_kind, basic_unit
675 : INTEGER, INTENT(in), OPTIONAL :: power
676 : REAL(kind=dp) :: res
677 :
678 : INTEGER :: my_power
679 :
680 5208477 : my_power = 1
681 5208477 : IF (PRESENT(power)) my_power = power
682 5208477 : IF (basic_unit == cp_units_none .AND. basic_kind /= cp_ukind_undef) THEN
683 0 : IF (basic_kind /= cp_units_none) THEN
684 : CALL cp_abort(__LOCATION__, &
685 : "unit not yet fully specified, unit of kind "// &
686 0 : TRIM(cp_to_string(basic_unit)))
687 : END IF
688 : END IF
689 5243545 : SELECT CASE (basic_kind)
690 : CASE (cp_ukind_undef)
691 873155 : SELECT CASE (basic_unit)
692 : CASE (cp_units_none)
693 35068 : res = value
694 : CASE default
695 35068 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
696 : END SELECT
697 : CASE (cp_ukind_energy)
698 3317629 : SELECT CASE (basic_unit)
699 : CASE (cp_units_hartree, cp_units_au)
700 160236 : res = value
701 : CASE (cp_units_wavenum)
702 0 : res = wavenumbers**(-my_power)*value
703 : CASE (cp_units_joule)
704 0 : res = joule**(-my_power)*value
705 : CASE (cp_units_kcalmol)
706 214608 : res = kcalmol**(-my_power)*value
707 : CASE (cp_units_kjmol)
708 2248 : res = kjmol**(-my_power)*value
709 : CASE (cp_units_jmol)
710 0 : res = (kjmol*1.0E+3_dp)**(-my_power)*value
711 : CASE (cp_units_Ry)
712 23209 : res = 0.5_dp**my_power*value
713 : CASE (cp_units_eV)
714 414050 : res = evolt**(-my_power)*value
715 : CASE (cp_units_keV)
716 8 : res = (1.0E-3_dp*evolt)**(-my_power)*value
717 : CASE (cp_units_k)
718 23728 : res = kelvin**(-my_power)*value
719 : CASE default
720 838087 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
721 : END SELECT
722 : CASE (cp_ukind_length)
723 363871 : SELECT CASE (basic_unit)
724 : CASE (cp_units_bohr, cp_units_au)
725 174525 : res = value
726 : CASE (cp_units_m)
727 58 : res = value*(1.0E10_dp*bohr)**my_power
728 : CASE (cp_units_pm)
729 2 : res = value*(0.01_dp*bohr)**my_power
730 : CASE (cp_units_nm)
731 10075 : res = value*(10.0_dp*bohr)**my_power
732 : CASE (cp_units_angstrom)
733 2972733 : res = value*bohr**my_power
734 : CASE default
735 3157393 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
736 : END SELECT
737 : CASE (cp_ukind_temperature)
738 232013 : SELECT CASE (basic_unit)
739 : CASE (cp_units_k)
740 189346 : res = kelvin**(-my_power)*value
741 : CASE (cp_units_au)
742 0 : res = value
743 : CASE default
744 189346 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
745 : END SELECT
746 : CASE (cp_ukind_pressure)
747 251422 : SELECT CASE (basic_unit)
748 : CASE (cp_units_bar)
749 41001 : res = bar**(-my_power)*value
750 : CASE (cp_units_atm)
751 0 : res = atm**(-my_power)*value
752 : CASE (cp_units_kbar)
753 18 : res = (1.0E-3_dp*bar)**(-my_power)*value
754 : CASE (cp_units_Pa)
755 1358 : res = pascal**(-my_power)*value
756 : CASE (cp_units_MPa)
757 0 : res = (1.0E-6_dp*pascal)**(-my_power)*value
758 : CASE (cp_units_GPa)
759 290 : res = (1.0E-9_dp*pascal)**(-my_power)*value
760 : CASE (cp_units_au)
761 0 : res = value
762 : CASE default
763 42667 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
764 : END SELECT
765 : CASE (cp_ukind_angle)
766 805316 : SELECT CASE (basic_unit)
767 : CASE (cp_units_rad)
768 70127 : res = value
769 : CASE (cp_units_deg)
770 140294 : res = value*(radians)**my_power
771 : CASE default
772 210421 : CPABORT("unknown angle unit:"//TRIM(cp_to_string(basic_unit)))
773 : END SELECT
774 : CASE (cp_ukind_time)
775 215 : SELECT CASE (basic_unit)
776 : CASE (cp_units_s)
777 24 : res = value*seconds**(-my_power)
778 : CASE (cp_units_fs)
779 387182 : res = value*femtoseconds**(-my_power)
780 : CASE (cp_units_ps)
781 299032 : res = value*picoseconds**(-my_power)
782 : CASE (cp_units_au)
783 48917 : res = value
784 : CASE (cp_units_wn)
785 34 : res = (twopi*wavenumbers)**(my_power)/value
786 : CASE default
787 735189 : CPABORT("unknown time unit:"//TRIM(cp_to_string(basic_unit)))
788 : END SELECT
789 : CASE (cp_ukind_mass)
790 77 : SELECT CASE (basic_unit)
791 : CASE (cp_units_kg)
792 0 : res = e_mass**my_power*value
793 : CASE (cp_units_amu)
794 182 : res = massunit**my_power*value
795 : CASE (cp_units_m_e, cp_units_au)
796 9 : res = value
797 : CASE default
798 191 : CPABORT("unknown mass unit:"//TRIM(cp_to_string(basic_unit)))
799 : END SELECT
800 : CASE (cp_ukind_potential)
801 113 : SELECT CASE (basic_unit)
802 : CASE (cp_units_volt)
803 77 : res = evolt**(-my_power)*value
804 : CASE (cp_units_au)
805 0 : res = value
806 : CASE default
807 77 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
808 : END SELECT
809 : CASE (cp_ukind_force)
810 14 : SELECT CASE (basic_unit)
811 : CASE (cp_units_Newton)
812 12 : res = value*newton**(-my_power)
813 : CASE (cp_units_mNewton)
814 24 : res = value*(1.0E+3*newton)**(-my_power)
815 : CASE (cp_units_au)
816 0 : res = value
817 : CASE default
818 36 : CPABORT("unknown force unit:"//TRIM(cp_to_string(basic_unit)))
819 : END SELECT
820 : CASE (cp_ukind_efield)
821 2 : SELECT CASE (basic_unit)
822 : CASE (cp_units_volt_per_m)
823 2 : res = (1.0E+10_dp*evolt*bohr)**(-my_power)*value
824 : CASE (cp_units_volt_per_nm)
825 0 : res = (10.0_dp*evolt*bohr)**(-my_power)*value
826 : CASE (cp_units_volt_per_angstrom)
827 0 : res = (evolt*bohr)**(-my_power)*value
828 : CASE (cp_units_au)
829 0 : res = value
830 : CASE default
831 2 : CPABORT("unknown electric field unit:"//TRIM(cp_to_string(basic_unit)))
832 : END SELECT
833 : CASE (cp_ukind_none)
834 : CALL cp_abort(__LOCATION__, &
835 : "if the kind of the unit is none also unit must be undefined,not:" &
836 0 : //TRIM(cp_to_string(basic_unit)))
837 : CASE default
838 5208477 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
839 : END SELECT
840 5208477 : END FUNCTION cp_basic_unit_to_cp2k
841 :
842 : ! **************************************************************************************************
843 : !> \brief returns the label of the current basic unit
844 : !> \param basic_kind the kind of the unit of the value
845 : !> \param basic_unit the unit of the value
846 : !> \param power the power of the unit (defaults to 1)
847 : !> \param accept_undefined ...
848 : !> \return ...
849 : !> \author fawzi
850 : ! **************************************************************************************************
851 8354027 : FUNCTION cp_basic_unit_desc(basic_kind, basic_unit, power, accept_undefined) &
852 : RESULT(res)
853 : INTEGER, INTENT(in) :: basic_kind, basic_unit
854 : INTEGER, INTENT(in), OPTIONAL :: power
855 : LOGICAL, INTENT(in), OPTIONAL :: accept_undefined
856 : CHARACTER(len=cp_unit_basic_desc_length) :: res
857 :
858 : INTEGER :: a, my_power
859 : LOGICAL :: my_accept_undefined
860 :
861 8354027 : my_power = 1
862 8354027 : res = ""
863 8354027 : my_accept_undefined = .FALSE.
864 8354027 : IF (accept_undefined) my_accept_undefined = accept_undefined
865 8354027 : IF (PRESENT(power)) my_power = power
866 8354027 : IF (basic_unit == cp_units_none) THEN
867 311324 : IF (.NOT. my_accept_undefined .AND. basic_kind == cp_units_none) THEN
868 : CALL cp_abort(__LOCATION__, "unit not yet fully specified, unit of kind "// &
869 0 : TRIM(cp_to_string(basic_kind)))
870 : END IF
871 : END IF
872 8665351 : SELECT CASE (basic_kind)
873 : CASE (cp_ukind_undef)
874 1732497 : SELECT CASE (basic_unit)
875 : CASE (cp_units_none)
876 311324 : res = "internal_cp2k"
877 : CASE DEFAULT
878 : CALL cp_abort(__LOCATION__, &
879 : "unit not yet fully specified, unit of kind "// &
880 311324 : TRIM(res))
881 : END SELECT
882 : CASE (cp_ukind_energy)
883 4576917 : SELECT CASE (basic_unit)
884 : CASE (cp_units_hartree, cp_units_au)
885 330616 : res = "hartree"
886 : CASE (cp_units_wavenum)
887 0 : res = "wavenumber_e"
888 : CASE (cp_units_joule)
889 0 : res = "joule"
890 : CASE (cp_units_kcalmol)
891 234924 : res = "kcalmol"
892 : CASE (cp_units_kjmol)
893 2248 : res = "kjmol"
894 : CASE (cp_units_jmol)
895 0 : res = "jmol"
896 : CASE (cp_units_Ry)
897 51565 : res = "Ry"
898 : CASE (cp_units_eV)
899 681424 : res = "eV"
900 : CASE (cp_units_keV)
901 21700 : res = "keV"
902 : CASE (cp_units_k)
903 98696 : res = "K_e"
904 : CASE (cp_units_none)
905 0 : res = "energy"
906 0 : IF (.NOT. my_accept_undefined) THEN
907 : CALL cp_abort(__LOCATION__, &
908 : "unit not yet fully specified, unit of kind "// &
909 0 : TRIM(res))
910 : END IF
911 : CASE default
912 1421173 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
913 : END SELECT
914 : CASE (cp_ukind_length)
915 1092955 : SELECT CASE (basic_unit)
916 : CASE (cp_units_bohr, cp_units_au)
917 608133 : res = "bohr"
918 : CASE (cp_units_m)
919 24497 : res = "m"
920 : CASE (cp_units_pm)
921 2 : res = "pm"
922 : CASE (cp_units_nm)
923 11417 : res = "nm"
924 : CASE (cp_units_angstrom)
925 3602252 : res = "angstrom"
926 : CASE default
927 0 : res = "length"
928 4246301 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
929 : END SELECT
930 : CASE (cp_ukind_temperature)
931 572933 : SELECT CASE (basic_unit)
932 : CASE (cp_units_k)
933 484822 : res = "K"
934 : CASE (cp_units_au)
935 0 : res = "au_temp"
936 : CASE (cp_units_none)
937 0 : res = "temperature"
938 0 : IF (.NOT. my_accept_undefined) THEN
939 : CALL cp_abort(__LOCATION__, &
940 : "unit not yet fully specified, unit of kind "// &
941 0 : TRIM(res))
942 : END IF
943 : CASE default
944 484822 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
945 : END SELECT
946 : CASE (cp_ukind_pressure)
947 322536 : SELECT CASE (basic_unit)
948 : CASE (cp_units_bar)
949 83711 : res = "bar"
950 : CASE (cp_units_atm)
951 0 : res = "atm"
952 : CASE (cp_units_kbar)
953 18 : res = "kbar"
954 : CASE (cp_units_Pa)
955 2716 : res = "Pa"
956 : CASE (cp_units_MPa)
957 0 : res = "MPa"
958 : CASE (cp_units_GPa)
959 1666 : res = "GPa"
960 : CASE (cp_units_au)
961 0 : res = "au_p"
962 : CASE (cp_units_none)
963 0 : res = "pressure"
964 0 : IF (.NOT. my_accept_undefined) THEN
965 : CALL cp_abort(__LOCATION__, &
966 : "unit not yet fully specified, unit of kind "// &
967 0 : TRIM(res))
968 : END IF
969 : CASE default
970 88111 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
971 : END SELECT
972 : CASE (cp_ukind_angle)
973 1604237 : SELECT CASE (basic_unit)
974 : CASE (cp_units_rad)
975 87249 : res = "rad"
976 : CASE (cp_units_deg)
977 151576 : res = "deg"
978 : CASE (cp_units_none)
979 0 : res = "angle"
980 0 : IF (.NOT. my_accept_undefined) THEN
981 : CALL cp_abort(__LOCATION__, &
982 : "unit not yet fully specified, unit of kind "// &
983 0 : TRIM(res))
984 : END IF
985 : CASE default
986 238825 : CPABORT("unknown angle unit:"//TRIM(cp_to_string(basic_unit)))
987 : END SELECT
988 : CASE (cp_ukind_time)
989 61281 : SELECT CASE (basic_unit)
990 : CASE (cp_units_s)
991 38071 : res = "s"
992 : CASE (cp_units_fs)
993 1021278 : res = "fs"
994 : CASE (cp_units_ps)
995 365474 : res = "ps"
996 : CASE (cp_units_au)
997 92131 : res = "au_t"
998 : CASE (cp_units_wn)
999 34 : res = "wavenumber_t"
1000 : CASE (cp_units_none)
1001 0 : res = "time"
1002 0 : IF (.NOT. my_accept_undefined) THEN
1003 : CALL cp_abort(__LOCATION__, &
1004 : "unit not yet fully specified, unit of kind "// &
1005 0 : TRIM(res))
1006 : END IF
1007 : CASE default
1008 1516988 : CPABORT("unknown time unit:"//TRIM(cp_to_string(basic_unit)))
1009 : END SELECT
1010 : CASE (cp_ukind_mass)
1011 17750 : SELECT CASE (basic_unit)
1012 : CASE (cp_units_kg)
1013 0 : res = "kg"
1014 : CASE (cp_units_amu)
1015 1514 : res = "amu"
1016 : CASE (cp_units_m_e, cp_units_au)
1017 21696 : res = "m_e"
1018 : CASE (cp_units_none)
1019 0 : res = "mass"
1020 0 : IF (.NOT. my_accept_undefined) THEN
1021 : CALL cp_abort(__LOCATION__, &
1022 : "unit not yet fully specified, unit of kind "// &
1023 0 : TRIM(res))
1024 : END IF
1025 : CASE default
1026 23210 : CPABORT("unknown mass unit:"//TRIM(cp_to_string(basic_unit)))
1027 : END SELECT
1028 : CASE (cp_ukind_potential)
1029 20538 : SELECT CASE (basic_unit)
1030 : CASE (cp_units_volt)
1031 17750 : res = "volt"
1032 : CASE (cp_units_au)
1033 0 : res = "au_pot"
1034 : CASE (cp_units_none)
1035 0 : res = "potential"
1036 0 : IF (.NOT. my_accept_undefined) THEN
1037 : CALL cp_abort(__LOCATION__, &
1038 : "unit not yet fully specified, unit of kind "// &
1039 0 : TRIM(res))
1040 : END IF
1041 : CASE default
1042 17750 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
1043 : END SELECT
1044 : CASE (cp_ukind_force)
1045 2747 : SELECT CASE (basic_unit)
1046 : CASE (cp_units_Newton)
1047 12 : res = "N"
1048 : CASE (cp_units_mNewton)
1049 2776 : res = "mN"
1050 : CASE (cp_units_au)
1051 0 : res = "au_f"
1052 : CASE (cp_units_none)
1053 0 : res = "force"
1054 0 : IF (.NOT. my_accept_undefined) THEN
1055 : CALL cp_abort(__LOCATION__, &
1056 : "unit not yet fully specified, unit of kind "// &
1057 0 : TRIM(res))
1058 : END IF
1059 : CASE default
1060 2788 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
1061 : END SELECT
1062 : CASE (cp_ukind_efield)
1063 2735 : SELECT CASE (basic_unit)
1064 : CASE (cp_units_volt_per_m)
1065 2735 : res = "Vm-1"
1066 : CASE (cp_units_volt_per_nm)
1067 0 : res = "Vnm-1"
1068 : CASE (cp_units_volt_per_angstrom)
1069 0 : res = "Vangstrom-1"
1070 : CASE (cp_units_au)
1071 0 : res = "au_efield"
1072 : CASE (cp_units_none)
1073 0 : res = "electric field"
1074 0 : IF (.NOT. my_accept_undefined) THEN
1075 : CALL cp_abort(__LOCATION__, &
1076 : "unit not yet fully specified, unit of kind "// &
1077 0 : TRIM(res))
1078 : END IF
1079 : CASE default
1080 2735 : CPABORT("unknown efield unit:"//TRIM(cp_to_string(basic_unit)))
1081 : END SELECT
1082 : CASE (cp_ukind_none)
1083 : CALL cp_abort(__LOCATION__, &
1084 : "if the kind of the unit is none also unit must be undefined,not:" &
1085 0 : //TRIM(cp_to_string(basic_unit)))
1086 : CASE default
1087 8354027 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
1088 : END SELECT
1089 8354027 : IF (my_power /= 1) THEN
1090 1282359 : a = LEN_TRIM(res)
1091 1282359 : CPASSERT(LEN(res) - a >= 3)
1092 1282359 : WRITE (res(a + 1:), "('^',i3)") my_power
1093 1282359 : CALL compress(res, .TRUE.)
1094 : END IF
1095 8354027 : END FUNCTION cp_basic_unit_desc
1096 :
1097 : ! **************************************************************************************************
1098 : !> \brief returns the "name" of the given unit
1099 : !> \param unit the unit to describe
1100 : !> \param defaults defaults for the undefined units, optional
1101 : !> \param accept_undefined if defaults is not present or is not associated
1102 : !> whether undefined units should be accepted (defaults to false)
1103 : !> \return ...
1104 : !> \author fawzi
1105 : ! **************************************************************************************************
1106 7731594 : FUNCTION cp_unit_desc(unit, defaults, accept_undefined) &
1107 : RESULT(res)
1108 : TYPE(cp_unit_type), INTENT(IN) :: unit
1109 : TYPE(cp_unit_set_type), INTENT(IN), OPTIONAL :: defaults
1110 : LOGICAL, INTENT(in), OPTIONAL :: accept_undefined
1111 : CHARACTER(len=cp_unit_desc_length) :: res
1112 :
1113 : INTEGER :: i, my_unit, pos
1114 : LOGICAL :: check, has_defaults, my_accept_undefined
1115 :
1116 7731594 : res = ""
1117 7731594 : pos = 1
1118 7731594 : my_accept_undefined = .FALSE.
1119 7731594 : IF (PRESENT(accept_undefined)) my_accept_undefined = accept_undefined
1120 16085621 : DO i = 1, unit%n_kinds
1121 8354027 : CPASSERT(unit%kind_id(i) /= 0)
1122 8354027 : CPASSERT(pos < LEN(res))
1123 8354027 : my_unit = unit%unit_id(i)
1124 8354027 : has_defaults = .FALSE.
1125 8354027 : IF (PRESENT(defaults)) has_defaults = ASSOCIATED(defaults%units(1)%unit)
1126 8354027 : IF (my_unit == 0) THEN
1127 0 : IF (has_defaults) THEN
1128 0 : my_unit = defaults%units(unit%kind_id(i))%unit%unit_id(1)
1129 : ELSE
1130 0 : check = my_accept_undefined .OR. unit%kind_id(i) /= 0
1131 0 : CPASSERT(check)
1132 : END IF
1133 : END IF
1134 8354027 : IF (i > 1) THEN
1135 630581 : res(pos:pos) = "*"
1136 630581 : pos = pos + 1
1137 : END IF
1138 : res(pos:) = TRIM(cp_basic_unit_desc(basic_kind=unit%kind_id(i), &
1139 : basic_unit=my_unit, accept_undefined=my_accept_undefined, &
1140 8354027 : power=unit%power(i)))
1141 16085621 : pos = LEN_TRIM(res) + 1
1142 : END DO
1143 :
1144 7731594 : END FUNCTION cp_unit_desc
1145 :
1146 : ! **************************************************************************************************
1147 : !> \brief transform a value to the internal cp2k units
1148 : !> \param value the value to convert
1149 : !> \param unit the unit of the result
1150 : !> \param defaults the defaults unit for those that are left free
1151 : !> (cp_units_none)
1152 : !> \param power the power of the unit (defaults to 1)
1153 : !> \return ...
1154 : !> \author fawzi
1155 : ! **************************************************************************************************
1156 4672087 : FUNCTION cp_unit_to_cp2k1(value, unit, defaults, power) RESULT(res)
1157 : REAL(kind=dp), INTENT(in) :: value
1158 : TYPE(cp_unit_type), INTENT(IN) :: unit
1159 : TYPE(cp_unit_set_type), INTENT(IN), OPTIONAL :: defaults
1160 : INTEGER, INTENT(in), OPTIONAL :: power
1161 : REAL(kind=dp) :: res
1162 :
1163 : INTEGER :: i_unit, my_basic_unit, my_power
1164 :
1165 4672087 : my_power = 1
1166 4672087 : IF (PRESENT(power)) my_power = power
1167 4672087 : res = value
1168 9880564 : DO i_unit = 1, unit%n_kinds
1169 5208477 : CPASSERT(unit%kind_id(i_unit) > 0)
1170 5208477 : my_basic_unit = unit%unit_id(i_unit)
1171 5208477 : IF (my_basic_unit == 0 .AND. unit%kind_id(i_unit) /= cp_ukind_undef) THEN
1172 0 : CPASSERT(PRESENT(defaults))
1173 0 : CPASSERT(ASSOCIATED(defaults%units(unit%kind_id(i_unit))%unit))
1174 0 : my_basic_unit = defaults%units(unit%kind_id(i_unit))%unit%unit_id(1)
1175 : END IF
1176 : res = cp_basic_unit_to_cp2k(value=res, basic_unit=my_basic_unit, &
1177 : basic_kind=unit%kind_id(i_unit), &
1178 9880564 : power=my_power*unit%power(i_unit))
1179 : END DO
1180 4672087 : END FUNCTION cp_unit_to_cp2k1
1181 :
1182 : ! **************************************************************************************************
1183 : !> \brief converts from the internal cp2k units to the given unit
1184 : !> \param value the value to convert
1185 : !> \param unit the unit of the result
1186 : !> \param defaults the defaults unit for those that are left free
1187 : !> (cp_units_none)
1188 : !> \param power the power of the unit (defaults to 1)
1189 : !> \return ...
1190 : !> \author fawzi
1191 : ! **************************************************************************************************
1192 783418 : FUNCTION cp_unit_from_cp2k1(value, unit, defaults, power) RESULT(res)
1193 : REAL(kind=dp), INTENT(in) :: value
1194 : TYPE(cp_unit_type), INTENT(IN) :: unit
1195 : TYPE(cp_unit_set_type), INTENT(IN), OPTIONAL :: defaults
1196 : INTEGER, INTENT(in), OPTIONAL :: power
1197 : REAL(kind=dp) :: res
1198 :
1199 : INTEGER :: my_power
1200 :
1201 783418 : my_power = 1
1202 783418 : IF (PRESENT(power)) my_power = power
1203 783418 : IF (PRESENT(defaults)) THEN
1204 : res = cp_unit_to_cp2k1(value=value, unit=unit, defaults=defaults, &
1205 0 : power=-my_power)
1206 : ELSE
1207 783418 : res = cp_unit_to_cp2k1(value=value, unit=unit, power=-my_power)
1208 : END IF
1209 783418 : END FUNCTION cp_unit_from_cp2k1
1210 :
1211 : ! **************************************************************************************************
1212 : !> \brief converts to the internal cp2k units to the given unit
1213 : !> \param value the value to convert
1214 : !> \param unit_str the unit of the result as string
1215 : !> \param defaults the defaults unit for those that are left free
1216 : !> (cp_units_none)
1217 : !> \param power the power of the unit (defaults to 1)
1218 : !> \return ...
1219 : !> \author fawzi
1220 : ! **************************************************************************************************
1221 3772718 : FUNCTION cp_unit_to_cp2k(value, unit_str, defaults, power) RESULT(res)
1222 : REAL(kind=dp), INTENT(in) :: value
1223 : CHARACTER(len=*), INTENT(in) :: unit_str
1224 : TYPE(cp_unit_set_type), INTENT(IN), OPTIONAL :: defaults
1225 : INTEGER, INTENT(in), OPTIONAL :: power
1226 : REAL(kind=dp) :: res
1227 :
1228 : TYPE(cp_unit_type) :: my_unit
1229 :
1230 3772718 : CALL cp_unit_create(my_unit, unit_str)
1231 3772718 : IF (PRESENT(defaults)) THEN
1232 : res = cp_unit_to_cp2k1(value=value, unit=my_unit, defaults=defaults, &
1233 0 : power=power)
1234 : ELSE
1235 3772718 : res = cp_unit_to_cp2k1(value=value, unit=my_unit, power=power)
1236 : END IF
1237 3772718 : CALL cp_unit_release(my_unit)
1238 101863386 : END FUNCTION cp_unit_to_cp2k
1239 :
1240 : ! **************************************************************************************************
1241 : !> \brief converts from the internal cp2k units to the given unit
1242 : !> \param value the value to convert
1243 : !> \param unit_str the unit of the result as string
1244 : !> \param defaults the defaults unit for those that are left free
1245 : !> (cp_units_none)
1246 : !> \param power the power of the unit (defaults to 1)
1247 : !> \return ...
1248 : !> \author fawzi
1249 : ! **************************************************************************************************
1250 584086 : FUNCTION cp_unit_from_cp2k(value, unit_str, defaults, power) RESULT(res)
1251 : REAL(kind=dp), INTENT(in) :: value
1252 : CHARACTER(len=*), INTENT(in) :: unit_str
1253 : TYPE(cp_unit_set_type), INTENT(IN), OPTIONAL :: defaults
1254 : INTEGER, INTENT(in), OPTIONAL :: power
1255 : REAL(kind=dp) :: res
1256 :
1257 : TYPE(cp_unit_type) :: my_unit
1258 :
1259 584086 : CALL cp_unit_create(my_unit, unit_str)
1260 584086 : IF (PRESENT(defaults)) THEN
1261 : res = cp_unit_from_cp2k1(value=value, unit=my_unit, defaults=defaults, &
1262 0 : power=power)
1263 : ELSE
1264 584086 : res = cp_unit_from_cp2k1(value=value, unit=my_unit, power=power)
1265 : END IF
1266 584086 : CALL cp_unit_release(my_unit)
1267 15770322 : END FUNCTION cp_unit_from_cp2k
1268 :
1269 : ! **************************************************************************************************
1270 : !> \brief returs true if the two units are compatible
1271 : !> \param ref_unit ...
1272 : !> \param unit ...
1273 : !> \return ...
1274 : !> \author Teodoro Laino [tlaino] - 11.2007 - University of Zurich
1275 : ! **************************************************************************************************
1276 115951 : FUNCTION cp_unit_compatible(ref_unit, unit) RESULT(res)
1277 : TYPE(cp_unit_type), INTENT(IN) :: ref_unit, unit
1278 : LOGICAL :: res
1279 :
1280 : INTEGER :: i
1281 :
1282 115951 : res = .TRUE.
1283 1043559 : DO i = 1, SIZE(ref_unit%kind_id)
1284 927608 : IF (ref_unit%kind_id(i) == unit%kind_id(i)) CYCLE
1285 5248 : IF ((ref_unit%kind_id(1) == cp_ukind_undef) .AND. (ALL(ref_unit%kind_id(2:) == cp_ukind_none))) CYCLE
1286 : res = .FALSE.
1287 1043559 : EXIT
1288 : END DO
1289 :
1290 115951 : END FUNCTION cp_unit_compatible
1291 :
1292 : ! **************************************************************************************************
1293 : !> \brief initializes the given unit set
1294 : !> \param unit_set the set to initialize
1295 : !> \param name the name of the set, used for the dafault initialization of
1296 : !> the various units
1297 : !> \author fawzi
1298 : ! **************************************************************************************************
1299 141708 : SUBROUTINE cp_unit_set_create(unit_set, name)
1300 : TYPE(cp_unit_set_type), INTENT(OUT) :: unit_set
1301 : CHARACTER(len=*), INTENT(in) :: name
1302 :
1303 : CHARACTER(len=cp_unit_desc_length) :: my_name
1304 : INTEGER :: i
1305 :
1306 11809 : my_name = name
1307 11809 : CALL uppercase(my_name)
1308 :
1309 141708 : DO i = 1, cp_ukind_max
1310 129899 : NULLIFY (unit_set%units(i)%unit)
1311 3259284 : ALLOCATE (unit_set%units(i)%unit)
1312 : END DO
1313 141708 : DO i = 1, cp_ukind_max
1314 11809 : SELECT CASE (name)
1315 : CASE ('ATOM', 'ATOMIC', 'INTERNAL', 'CP2K')
1316 0 : IF (i == cp_ukind_angle) THEN
1317 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], &
1318 0 : unit_id=[cp_units_rad], power=[1])
1319 : ELSE
1320 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], &
1321 0 : unit_id=[cp_units_au], power=[1])
1322 : END IF
1323 : CASE ('OUTPUT')
1324 11809 : SELECT CASE (i)
1325 : CASE (cp_ukind_undef)
1326 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_none], &
1327 23618 : power=[1])
1328 : CASE (cp_ukind_energy)
1329 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_hartree], &
1330 23618 : power=[1])
1331 : CASE (cp_ukind_length)
1332 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_angstrom], &
1333 23618 : power=[1])
1334 : CASE (cp_ukind_temperature)
1335 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_k], &
1336 23618 : power=[1])
1337 : CASE (cp_ukind_angle)
1338 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_deg], &
1339 23618 : power=[1])
1340 : CASE (cp_ukind_pressure)
1341 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_bar], &
1342 23618 : power=[1])
1343 : CASE (cp_ukind_time)
1344 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_fs], &
1345 23618 : power=[1])
1346 : CASE (cp_ukind_mass)
1347 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_amu], &
1348 23618 : power=[1])
1349 : CASE (cp_ukind_potential)
1350 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_volt], &
1351 23618 : power=[1])
1352 : CASE (cp_ukind_force)
1353 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_newton], &
1354 23618 : power=[1])
1355 : CASE (cp_ukind_efield)
1356 : CALL cp_unit_create2(unit_set%units(i)%unit, kind_id=[i], unit_id=[cp_units_volt_per_m], &
1357 23618 : power=[1])
1358 : CASE default
1359 0 : CPABORT("unhandled unit type "//TRIM(cp_to_string(i)))
1360 129899 : EXIT
1361 : END SELECT
1362 : CASE default
1363 129899 : CPABORT('unknown parameter set name '//TRIM(name))
1364 : END SELECT
1365 : END DO
1366 11809 : END SUBROUTINE cp_unit_set_create
1367 :
1368 : ! **************************************************************************************************
1369 : !> \brief releases the given unit set
1370 : !> \param unit_set the unit set to release
1371 : !> \author fawzi
1372 : ! **************************************************************************************************
1373 11809 : SUBROUTINE cp_unit_set_release(unit_set)
1374 : TYPE(cp_unit_set_type), INTENT(INOUT) :: unit_set
1375 :
1376 : INTEGER :: i
1377 :
1378 141708 : DO i = 1, SIZE(unit_set%units)
1379 129899 : CALL cp_unit_release(unit_set%units(i)%unit)
1380 141708 : DEALLOCATE (unit_set%units(i)%unit)
1381 : END DO
1382 :
1383 11809 : END SUBROUTINE cp_unit_set_release
1384 :
1385 : ! **************************************************************************************************
1386 : !> \brief Exports all available units as XML.
1387 : !> \param iw ...
1388 : !> \author Ole Schuett
1389 : ! **************************************************************************************************
1390 0 : SUBROUTINE export_units_as_xml(iw)
1391 : INTEGER, INTENT(IN) :: iw
1392 :
1393 : CALL format_units_as_xml("energy", s2a("hartree", "wavenumber_e", "joule", "kcalmol", &
1394 0 : "kjmol", "Ry", "eV", "keV", "K_e"), iw)
1395 0 : CALL format_units_as_xml("length", s2a("bohr", "m", "pm", "nm", "angstrom"), iw)
1396 0 : CALL format_units_as_xml("temperature", s2a("K", "au_temp"), iw)
1397 0 : CALL format_units_as_xml("pressure", s2a("bar", "atm", "kbar", "Pa", "MPa", "GPa", "au_p"), iw)
1398 0 : CALL format_units_as_xml("angle", s2a("rad", "deg"), iw)
1399 0 : CALL format_units_as_xml("time", s2a("s", "fs", "ps", "au_t", "wavenumber_t"), iw)
1400 0 : CALL format_units_as_xml("mass", s2a("kg", "amu", "m_e"), iw)
1401 0 : CALL format_units_as_xml("potential", s2a("volt", "au_pot"), iw)
1402 0 : CALL format_units_as_xml("force", s2a("N", "Newton", "mN", "mNewton", "au_f"), iw)
1403 : CALL format_units_as_xml("efield", s2a("Vm-1", "Vnm-1", "VA-1", "volt_per_m", "volt_per_nm", &
1404 0 : "volt_per_angstrom", "au_efield"), iw)
1405 :
1406 0 : END SUBROUTINE export_units_as_xml
1407 :
1408 : ! **************************************************************************************************
1409 : !> \brief Format units as xml.
1410 : !> \param unit_kind ...
1411 : !> \param units_set ...
1412 : !> \param iw ...
1413 : !> \author Ole Schuett
1414 : ! **************************************************************************************************
1415 0 : SUBROUTINE format_units_as_xml(unit_kind, units_set, iw)
1416 : CHARACTER(LEN=*), INTENT(IN) :: unit_kind
1417 : CHARACTER(LEN=*), DIMENSION(:), INTENT(IN) :: units_set
1418 : INTEGER, INTENT(IN) :: iw
1419 :
1420 : INTEGER :: i
1421 :
1422 0 : WRITE (iw, FMT='(T2,A)') '<UNIT_KIND name="'//TRIM(unit_kind)//'">'
1423 0 : DO i = 1, SIZE(units_set)
1424 0 : WRITE (iw, FMT='(T3,A)') '<UNIT>'//TRIM(units_set(i))//'</UNIT>'
1425 : END DO
1426 0 : WRITE (iw, FMT='(T3,A)') '<UNIT>'//TRIM(unit_kind)//'</UNIT>' ! internal unit
1427 0 : WRITE (iw, FMT='(T2,A)') '</UNIT_KIND>'
1428 0 : END SUBROUTINE format_units_as_xml
1429 :
1430 0 : END MODULE cp_units
|
