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 563745150 : 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 22549806 : 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 202948254 : unit_id = cp_units_none
179 22549806 : kind_id = cp_ukind_none
180 22549806 : power = 0
181 22549806 : i_low = 1
182 22549806 : i_high = 1
183 22549806 : len_string = LEN(string)
184 22549806 : i_unit = 0
185 22549806 : next_power = 1
186 22549806 : DO WHILE (i_low < len_string)
187 21502017 : IF (string(i_low:i_low) /= ' ') EXIT
188 22549806 : i_low = i_low + 1
189 : END DO
190 : i_high = i_low
191 147037965 : DO WHILE (i_high <= len_string)
192 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
193 127095892 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
194 144430287 : i_high = i_high + 1
195 : END DO
196 : DO
197 23527276 : IF (i_high <= i_low .OR. i_low > len_string) EXIT
198 23464343 : i_unit = i_unit + 1
199 23464343 : 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 23464343 : unit_string = string(i_low:i_high - 1)
205 23464343 : CALL uppercase(unit_string)
206 24496567 : SELECT CASE (TRIM(unit_string))
207 : CASE ("INTERNAL_CP2K")
208 1032224 : unit_id(i_unit) = cp_units_none
209 1032224 : kind_id(i_unit) = cp_ukind_undef
210 : CASE ("HARTREE")
211 1084704 : unit_id(i_unit) = cp_units_hartree
212 1084704 : kind_id(i_unit) = cp_ukind_energy
213 : CASE ("AU_E")
214 172870 : unit_id(i_unit) = cp_units_au
215 172870 : 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 390100 : unit_id(i_unit) = cp_units_kcalmol
224 390100 : 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 334050 : unit_id(i_unit) = cp_units_Ry
233 334050 : kind_id(i_unit) = cp_ukind_energy
234 : CASE ("EV")
235 3355994 : unit_id(i_unit) = cp_units_eV
236 3355994 : kind_id(i_unit) = cp_ukind_energy
237 : CASE ("KEV")
238 30809 : unit_id(i_unit) = cp_units_keV
239 30809 : kind_id(i_unit) = cp_ukind_energy
240 : CASE ("K_E")
241 345095 : unit_id(i_unit) = cp_units_k
242 345095 : 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 3045558 : unit_id(i_unit) = cp_units_bohr
251 3045558 : kind_id(i_unit) = cp_ukind_length
252 : CASE ("M")
253 51856 : unit_id(i_unit) = cp_units_m
254 51856 : 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 20540 : unit_id(i_unit) = cp_units_nm
260 20540 : kind_id(i_unit) = cp_ukind_length
261 : CASE ("ANGSTROM")
262 9328239 : unit_id(i_unit) = cp_units_angstrom
263 9328239 : 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 1020689 : unit_id(i_unit) = cp_units_k
269 1020689 : 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 156663 : unit_id(i_unit) = cp_units_bar
281 156663 : 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 20972 : unit_id(i_unit) = cp_units_Pa
287 20972 : 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 10791 : unit_id(i_unit) = cp_units_GPa
293 10791 : 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 251541 : unit_id(i_unit) = cp_units_rad
302 251541 : kind_id(i_unit) = cp_ukind_angle
303 : CASE ("DEG")
304 306534 : unit_id(i_unit) = cp_units_deg
305 306534 : 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 156704 : unit_id(i_unit) = cp_units_s
311 156704 : kind_id(i_unit) = cp_ukind_time
312 : CASE ("FS")
313 1577116 : unit_id(i_unit) = cp_units_fs
314 1577116 : kind_id(i_unit) = cp_ukind_time
315 : CASE ("PS")
316 438384 : unit_id(i_unit) = cp_units_ps
317 438384 : 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 110459 : unit_id(i_unit) = cp_units_au
323 110459 : 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 10638 : unit_id(i_unit) = cp_units_amu
332 10638 : 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 30805 : unit_id(i_unit) = cp_units_au
338 30805 : 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 136402 : unit_id(i_unit) = cp_units_volt
344 136402 : 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 21026 : unit_id(i_unit) = cp_units_mNewton
356 21026 : 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 20988 : unit_id(i_unit) = cp_units_volt_per_m
365 20988 : 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 23464343 : CPABORT("Unknown unit: "//string(i_low:i_high - 1))
385 : END SELECT
386 23464343 : power(i_unit) = next_power
387 : ! parse op
388 23464343 : i_low = i_high
389 23487357 : DO WHILE (i_low <= len_string)
390 3449308 : IF (string(i_low:i_low) /= ' ') EXIT
391 23487357 : i_low = i_low + 1
392 : END DO
393 : i_high = i_low
394 23464343 : DO WHILE (i_high <= len_string)
395 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
396 3426294 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
397 0 : i_high = i_high + 1
398 : END DO
399 23464343 : IF (i_high < i_low .OR. i_low > len_string) EXIT
400 :
401 3426294 : IF (i_high <= len_string) THEN
402 3426294 : IF (string(i_low:i_high) == '^') THEN
403 2491005 : i_low = i_high + 1
404 2491005 : DO WHILE (i_low <= len_string)
405 2491005 : IF (string(i_low:i_low) /= ' ') EXIT
406 2491005 : i_low = i_low + 1
407 : END DO
408 : i_high = i_low
409 7292981 : DO WHILE (i_high <= len_string)
410 2448824 : SELECT CASE (string(i_high:i_high))
411 : CASE ('+', '-', '0', '1', '2', '3', '4', '5', '6', '7', '8', '9')
412 4801976 : i_high = i_high + 1
413 : CASE default
414 4844157 : EXIT
415 : END SELECT
416 : END DO
417 2491005 : 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 2491005 : formatstr = "(i"//cp_to_string(i_high - i_low + 1)//")"
422 : READ (string(i_low:i_high - 1), formatstr) &
423 2491005 : next_power
424 2491005 : power(i_unit) = power(i_unit)*next_power
425 : ! next op
426 2491005 : i_low = i_high
427 2495925 : DO WHILE (i_low < len_string)
428 47084 : IF (string(i_low:i_low) /= ' ') EXIT
429 2495925 : i_low = i_low + 1
430 : END DO
431 : i_high = i_low
432 2492199 : DO WHILE (i_high <= len_string)
433 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
434 43173 : 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 3426294 : IF (i_low > len_string) EXIT
440 977470 : next_power = 1
441 23527276 : IF (i_high <= len_string) THEN
442 977268 : IF (string(i_low:i_high) == "*" .OR. string(i_low:i_high) == '/') THEN
443 977251 : IF (string(i_low:i_high) == '/') next_power = -1
444 977251 : i_low = i_high + 1
445 977251 : DO WHILE (i_low <= len_string)
446 977251 : IF (string(i_low:i_low) /= ' ') EXIT
447 977251 : i_low = i_low + 1
448 : END DO
449 : i_high = i_low
450 5218232 : DO WHILE (i_high <= len_string)
451 : IF (string(i_high:i_high) == ' ' .OR. string(i_high:i_high) == '^' .OR. &
452 5059599 : string(i_high:i_high) == '*' .OR. string(i_high:i_high) == '/') EXIT
453 4399633 : 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 22549806 : power=power)
460 22549806 : desc = cp_unit_desc(unit)
461 22549806 : 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 634901148 : 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 22675041 : CPASSERT(SIZE(kind_id) <= cp_unit_max_kinds)
482 22675041 : CPASSERT(SIZE(unit_id) <= cp_unit_max_kinds)
483 203198724 : unit%kind_id(1:SIZE(kind_id)) = kind_id
484 23551686 : unit%kind_id(SIZE(kind_id) + 1:) = cp_ukind_none
485 203198724 : unit%unit_id(1:SIZE(unit_id)) = unit_id
486 46226727 : unit%unit_id(SIZE(unit_id):) = cp_units_none
487 22675041 : IF (PRESENT(power)) THEN
488 203198724 : unit%power(1:SIZE(power)) = power
489 23551686 : unit%power(SIZE(power) + 1:) = 0
490 204075369 : DO i = 1, SIZE(unit%power)
491 204075369 : IF (unit%power(i) == 0) THEN
492 157810750 : unit%kind_id(i) = cp_ukind_none
493 157810750 : 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 22675041 : unit%n_kinds = 0
509 204075369 : DO i = 1, SIZE(unit%kind_id)
510 : ! find max and compress in the rest
511 : DO
512 181400328 : max_kind = unit%kind_id(i)
513 181400328 : max_pos = i
514 181400328 : repeat = .FALSE.
515 816301476 : DO j = i + 1, SIZE(unit%kind_id)
516 816301476 : IF (unit%kind_id(j) >= max_kind) THEN
517 471707748 : 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 471707748 : 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 935217 : max_kind = unit%kind_id(j)
534 935217 : max_pos = j
535 : END IF
536 : END IF
537 : END DO
538 181400328 : IF (.NOT. repeat) EXIT
539 : END DO
540 181400328 : IF (max_kind /= 0) unit%n_kinds = unit%n_kinds + 1
541 : ! put the max at pos i
542 181400328 : IF (max_pos /= i) THEN
543 914241 : unit%kind_id(max_pos) = unit%kind_id(i)
544 914241 : unit%kind_id(i) = max_kind
545 914241 : max_kind = unit%unit_id(max_pos)
546 914241 : unit%unit_id(max_pos) = unit%unit_id(i)
547 914241 : unit%unit_id(i) = max_kind
548 914241 : max_kind = unit%power(max_pos)
549 914241 : unit%power(max_pos) = unit%power(i)
550 914241 : unit%power(i) = max_kind
551 : END IF
552 : ! check unit
553 : CALL cp_basic_unit_check(basic_kind=unit%kind_id(i), &
554 204075369 : basic_unit=unit%unit_id(i))
555 : END DO
556 22675041 : 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 22675041 : ELEMENTAL SUBROUTINE cp_unit_release(unit)
566 : TYPE(cp_unit_type), INTENT(IN) :: unit
567 :
568 : MARK_USED(unit)
569 :
570 22675041 : 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 181400328 : SUBROUTINE cp_basic_unit_check(basic_kind, basic_unit)
579 : INTEGER, INTENT(in) :: basic_kind, basic_unit
580 :
581 182443937 : SELECT CASE (basic_kind)
582 : CASE (cp_ukind_undef)
583 6770864 : SELECT CASE (basic_unit)
584 : CASE (cp_units_none)
585 : CASE default
586 1043609 : CPABORT("unknown undef unit:"//TRIM(cp_to_string(basic_unit)))
587 : END SELECT
588 : CASE (cp_ukind_energy)
589 18185113 : 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 5727255 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
595 : END SELECT
596 : CASE (cp_ukind_length)
597 13489932 : 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 12457858 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
602 : END SELECT
603 : CASE (cp_ukind_temperature)
604 1231903 : SELECT CASE (basic_unit)
605 : CASE (cp_units_k, cp_units_au, cp_units_none)
606 : CASE default
607 1032074 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
608 : END SELECT
609 : CASE (cp_ukind_pressure)
610 769289 : 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 199829 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
614 : END SELECT
615 : CASE (cp_ukind_angle)
616 2863542 : SELECT CASE (basic_unit)
617 : CASE (cp_units_rad, cp_units_deg, cp_units_none)
618 : CASE default
619 569460 : CPABORT("unknown angle unit:"//TRIM(cp_to_string(basic_unit)))
620 : END SELECT
621 : CASE (cp_ukind_time)
622 2346910 : 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 2294082 : CPABORT("unknown time unit:"//TRIM(cp_to_string(basic_unit)))
626 : END SELECT
627 : CASE (cp_ukind_mass)
628 200615 : 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 52828 : CPABORT("unknown mass unit:"//TRIM(cp_to_string(basic_unit)))
632 : END SELECT
633 : CASE (cp_ukind_potential)
634 180210 : SELECT CASE (basic_unit)
635 : CASE (cp_units_volt, cp_units_au, cp_units_none)
636 : CASE default
637 147787 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
638 : END SELECT
639 : CASE (cp_ukind_force)
640 64796 : SELECT CASE (basic_unit)
641 : CASE (cp_units_Newton, cp_units_mNewton, cp_units_au, cp_units_none)
642 : CASE default
643 32423 : CPABORT("unknown force unit:"//TRIM(cp_to_string(basic_unit)))
644 : END SELECT
645 : CASE (cp_ukind_efield)
646 157843123 : 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 32373 : CPABORT("unknown electric field unit:"//TRIM(cp_to_string(basic_unit)))
651 : END SELECT
652 : CASE (cp_ukind_none)
653 157810750 : 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 181400328 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
660 : END SELECT
661 181400328 : 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 8218492 : 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 8218492 : my_power = 1
681 8218492 : IF (PRESENT(power)) my_power = power
682 8218492 : 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 8290047 : SELECT CASE (basic_kind)
690 : CASE (cp_ukind_undef)
691 1910251 : SELECT CASE (basic_unit)
692 : CASE (cp_units_none)
693 71555 : res = value
694 : CASE default
695 71555 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
696 : END SELECT
697 : CASE (cp_ukind_energy)
698 4725728 : SELECT CASE (basic_unit)
699 : CASE (cp_units_hartree, cp_units_au)
700 178093 : 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 232864 : 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 49549 : res = 0.5_dp**my_power*value
713 : CASE (cp_units_eV)
714 1352270 : 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 23664 : res = kelvin**(-my_power)*value
719 : CASE default
720 1838696 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
721 : END SELECT
722 : CASE (cp_ukind_length)
723 865184 : SELECT CASE (basic_unit)
724 : CASE (cp_units_bohr, cp_units_au)
725 475734 : 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 10066 : res = value*(10.0_dp*bohr)**my_power
732 : CASE (cp_units_angstrom)
733 4061775 : res = value*bohr**my_power
734 : CASE default
735 4547635 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
736 : END SELECT
737 : CASE (cp_ukind_temperature)
738 477719 : SELECT CASE (basic_unit)
739 : CASE (cp_units_k)
740 389450 : res = kelvin**(-my_power)*value
741 : CASE (cp_units_au)
742 0 : res = value
743 : CASE default
744 389450 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
745 : END SELECT
746 : CASE (cp_ukind_pressure)
747 406232 : SELECT CASE (basic_unit)
748 : CASE (cp_units_bar)
749 77475 : 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 10486 : 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 88269 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
764 : END SELECT
765 : CASE (cp_ukind_angle)
766 1023947 : SELECT CASE (basic_unit)
767 : CASE (cp_units_rad)
768 70123 : res = value
769 : CASE (cp_units_deg)
770 258634 : res = value*(radians)**my_power
771 : CASE default
772 328757 : 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 569252 : res = value*femtoseconds**(-my_power)
780 : CASE (cp_units_ps)
781 335487 : res = value*picoseconds**(-my_power)
782 : CASE (cp_units_au)
783 49027 : res = value
784 : CASE (cp_units_wn)
785 34 : res = (twopi*wavenumbers)**(my_power)/value
786 : CASE default
787 953824 : 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 8218492 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
839 : END SELECT
840 8218492 : 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 23464343 : 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 23464343 : my_power = 1
862 23464343 : res = ""
863 23464343 : my_accept_undefined = .FALSE.
864 23464343 : IF (accept_undefined) my_accept_undefined = accept_undefined
865 23464343 : IF (PRESENT(power)) my_power = power
866 23464343 : IF (basic_unit == cp_units_none) THEN
867 1032224 : 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 24496567 : SELECT CASE (basic_kind)
873 : CASE (cp_ukind_undef)
874 6748094 : SELECT CASE (basic_unit)
875 : CASE (cp_units_none)
876 1032224 : res = "internal_cp2k"
877 : CASE DEFAULT
878 : CALL cp_abort(__LOCATION__, &
879 : "unit not yet fully specified, unit of kind "// &
880 1032224 : TRIM(res))
881 : END SELECT
882 : CASE (cp_ukind_energy)
883 13704047 : SELECT CASE (basic_unit)
884 : CASE (cp_units_hartree, cp_units_au)
885 1257574 : 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 390100 : 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 334050 : res = "Ry"
898 : CASE (cp_units_eV)
899 3355994 : res = "eV"
900 : CASE (cp_units_keV)
901 30809 : res = "keV"
902 : CASE (cp_units_k)
903 345095 : 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 5715870 : CPABORT("unknown energy unit:"//TRIM(cp_to_string(basic_unit)))
913 : END SELECT
914 : CASE (cp_ukind_length)
915 4066525 : SELECT CASE (basic_unit)
916 : CASE (cp_units_bohr, cp_units_au)
917 3045836 : res = "bohr"
918 : CASE (cp_units_m)
919 51856 : res = "m"
920 : CASE (cp_units_pm)
921 2 : res = "pm"
922 : CASE (cp_units_nm)
923 20540 : res = "nm"
924 : CASE (cp_units_angstrom)
925 9328239 : res = "angstrom"
926 : CASE default
927 0 : res = "length"
928 12446473 : CPABORT("unknown length unit:"//TRIM(cp_to_string(basic_unit)))
929 : END SELECT
930 : CASE (cp_ukind_temperature)
931 1209133 : SELECT CASE (basic_unit)
932 : CASE (cp_units_k)
933 1020689 : 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 1020689 : CPABORT("unknown temperature unit:"//TRIM(cp_to_string(basic_unit)))
945 : END SELECT
946 : CASE (cp_ukind_pressure)
947 714738 : SELECT CASE (basic_unit)
948 : CASE (cp_units_bar)
949 156663 : 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 20972 : res = "Pa"
956 : CASE (cp_units_MPa)
957 0 : res = "MPa"
958 : CASE (cp_units_GPa)
959 10791 : 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 188444 : CPABORT("unknown pressure unit:"//TRIM(cp_to_string(basic_unit)))
971 : END SELECT
972 : CASE (cp_ukind_angle)
973 2534238 : SELECT CASE (basic_unit)
974 : CASE (cp_units_rad)
975 251541 : res = "rad"
976 : CASE (cp_units_deg)
977 306534 : 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 558075 : CPABORT("unknown angle unit:"//TRIM(cp_to_string(basic_unit)))
987 : END SELECT
988 : CASE (cp_ukind_time)
989 198147 : SELECT CASE (basic_unit)
990 : CASE (cp_units_s)
991 156704 : res = "s"
992 : CASE (cp_units_fs)
993 1577116 : res = "fs"
994 : CASE (cp_units_ps)
995 438384 : res = "ps"
996 : CASE (cp_units_au)
997 110459 : 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 2282697 : CPABORT("unknown time unit:"//TRIM(cp_to_string(basic_unit)))
1009 : END SELECT
1010 : CASE (cp_ukind_mass)
1011 136402 : SELECT CASE (basic_unit)
1012 : CASE (cp_units_kg)
1013 0 : res = "kg"
1014 : CASE (cp_units_amu)
1015 10638 : res = "amu"
1016 : CASE (cp_units_m_e, cp_units_au)
1017 30805 : 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 41443 : CPABORT("unknown mass unit:"//TRIM(cp_to_string(basic_unit)))
1027 : END SELECT
1028 : CASE (cp_ukind_potential)
1029 157440 : SELECT CASE (basic_unit)
1030 : CASE (cp_units_volt)
1031 136402 : 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 136402 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
1043 : END SELECT
1044 : CASE (cp_ukind_force)
1045 21000 : SELECT CASE (basic_unit)
1046 : CASE (cp_units_Newton)
1047 12 : res = "N"
1048 : CASE (cp_units_mNewton)
1049 21026 : 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 21038 : CPABORT("unknown potential unit:"//TRIM(cp_to_string(basic_unit)))
1061 : END SELECT
1062 : CASE (cp_ukind_efield)
1063 20988 : SELECT CASE (basic_unit)
1064 : CASE (cp_units_volt_per_m)
1065 20988 : 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 20988 : 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 23464343 : CPABORT("unknown kind of unit:"//TRIM(cp_to_string(basic_kind)))
1088 : END SELECT
1089 23464343 : IF (my_power /= 1) THEN
1090 2628654 : a = LEN_TRIM(res)
1091 2628654 : CPASSERT(LEN(res) - a >= 3)
1092 2628654 : WRITE (res(a + 1:), "('^',i3)") my_power
1093 2628654 : CALL compress(res, .TRUE.)
1094 : END IF
1095 23464343 : 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 22549806 : 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 22549806 : res = ""
1117 22549806 : pos = 1
1118 22549806 : my_accept_undefined = .FALSE.
1119 22549806 : IF (PRESENT(accept_undefined)) my_accept_undefined = accept_undefined
1120 46014149 : DO i = 1, unit%n_kinds
1121 23464343 : CPASSERT(unit%kind_id(i) /= 0)
1122 23464343 : CPASSERT(pos < LEN(res))
1123 23464343 : my_unit = unit%unit_id(i)
1124 23464343 : has_defaults = .FALSE.
1125 23464343 : IF (PRESENT(defaults)) has_defaults = ASSOCIATED(defaults%units(1)%unit)
1126 23464343 : 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 23464343 : IF (i > 1) THEN
1135 977453 : res(pos:pos) = "*"
1136 977453 : 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 23464343 : power=unit%power(i)))
1141 46014149 : pos = LEN_TRIM(res) + 1
1142 : END DO
1143 :
1144 22549806 : 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 7672968 : 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 7672968 : my_power = 1
1166 7672968 : IF (PRESENT(power)) my_power = power
1167 7672968 : res = value
1168 15891460 : DO i_unit = 1, unit%n_kinds
1169 8218492 : CPASSERT(unit%kind_id(i_unit) > 0)
1170 8218492 : my_basic_unit = unit%unit_id(i_unit)
1171 8218492 : 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 15891460 : power=my_power*unit%power(i_unit))
1179 : END DO
1180 7672968 : 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 776342 : 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 776342 : my_power = 1
1202 776342 : IF (PRESENT(power)) my_power = power
1203 776342 : IF (PRESENT(defaults)) THEN
1204 : res = cp_unit_to_cp2k1(value=value, unit=unit, defaults=defaults, &
1205 0 : power=-my_power)
1206 : ELSE
1207 776342 : res = cp_unit_to_cp2k1(value=value, unit=unit, power=-my_power)
1208 : END IF
1209 776342 : 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 6783303 : 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 6783303 : CALL cp_unit_create(my_unit, unit_str)
1231 6783303 : IF (PRESENT(defaults)) THEN
1232 : res = cp_unit_to_cp2k1(value=value, unit=my_unit, defaults=defaults, &
1233 0 : power=power)
1234 : ELSE
1235 6783303 : res = cp_unit_to_cp2k1(value=value, unit=my_unit, power=power)
1236 : END IF
1237 6783303 : CALL cp_unit_release(my_unit)
1238 183149181 : 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 577420 : 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 577420 : CALL cp_unit_create(my_unit, unit_str)
1260 577420 : IF (PRESENT(defaults)) THEN
1261 : res = cp_unit_from_cp2k1(value=value, unit=my_unit, defaults=defaults, &
1262 0 : power=power)
1263 : ELSE
1264 577420 : res = cp_unit_from_cp2k1(value=value, unit=my_unit, power=power)
1265 : END IF
1266 577420 : CALL cp_unit_release(my_unit)
1267 15590340 : 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 113323 : 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 113323 : res = .TRUE.
1283 1019907 : DO i = 1, SIZE(ref_unit%kind_id)
1284 906584 : IF (ref_unit%kind_id(i) == unit%kind_id(i)) CYCLE
1285 5200 : IF ((ref_unit%kind_id(1) == cp_ukind_undef) .AND. (ALL(ref_unit%kind_id(2:) == cp_ukind_none))) CYCLE
1286 : res = .FALSE.
1287 1019907 : EXIT
1288 : END DO
1289 :
1290 113323 : 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 136620 : 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 11385 : my_name = name
1307 11385 : CALL uppercase(my_name)
1308 :
1309 136620 : DO i = 1, cp_ukind_max
1310 125235 : NULLIFY (unit_set%units(i)%unit)
1311 3142260 : ALLOCATE (unit_set%units(i)%unit)
1312 : END DO
1313 136620 : DO i = 1, cp_ukind_max
1314 11385 : 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 11385 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : 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 22770 : power=[1])
1358 : CASE default
1359 0 : CPABORT("unhandled unit type "//TRIM(cp_to_string(i)))
1360 125235 : EXIT
1361 : END SELECT
1362 : CASE default
1363 125235 : CPABORT('unknown parameter set name '//TRIM(name))
1364 : END SELECT
1365 : END DO
1366 11385 : 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 11385 : SUBROUTINE cp_unit_set_release(unit_set)
1374 : TYPE(cp_unit_set_type), INTENT(INOUT) :: unit_set
1375 :
1376 : INTEGER :: i
1377 :
1378 136620 : DO i = 1, SIZE(unit_set%units)
1379 125235 : CALL cp_unit_release(unit_set%units(i)%unit)
1380 136620 : DEALLOCATE (unit_set%units(i)%unit)
1381 : END DO
1382 :
1383 11385 : 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
|
