20 integer,
parameter :: max_units = 200
22 integer,
parameter :: unit_offset = 10
24 logical,
save :: unit_used(max_units) = .false.
27 integer,
parameter :: PMC_UTIL_CONVERT_STRING_LEN = 100
29 integer,
parameter :: PMC_MAX_FILENAME_LEN = 300
36 subroutine warn_msg(code, warning_msg, already_warned)
39 integer,
intent(in) :: code
41 character(len=*),
intent(in) :: warning_msg
43 logical,
intent(inout),
optional :: already_warned
45 if (present(already_warned))
then
46 if (already_warned)
return
48 already_warned = .true.
51 //
'): ' // trim(warning_msg)
61 integer,
intent(in) :: code
63 logical,
intent(in) :: condition_ok
65 character(len=*),
intent(in) :: warning_msg
67 if (.not. condition_ok)
then
79 integer,
intent(in) :: code
81 logical,
intent(in) :: condition_ok
83 character(len=*),
intent(in) :: error_msg
87 if (.not. condition_ok)
then
89 //
'): ' // trim(error_msg)
91 call mpi_abort(mpi_comm_world, code, ierr)
105 integer,
intent(in) :: code
107 logical,
intent(in) :: condition_ok
109 if (.not. condition_ok)
then
113 call
assert_msg(code, condition_ok,
'assertion failed')
124 integer,
intent(in) :: code
126 call
assert(code, .false.)
136 integer,
intent(in) :: code
138 character(len=*),
intent(in) :: error_msg
154 if (.not. unit_used(i))
then
159 if (.not. found_unit)
then
161 'no more units available - need to free_unit()')
163 unit_used(i) = .true.
173 integer,
intent(in) :: unit
175 unit_used(unit - unit_offset) = .false.
187 character(len=*),
intent(in) :: filename
189 integer,
intent(out) :: unit
194 open(unit=unit, file=filename, status=
'old', action=
'read', iostat=ios)
195 call
assert_msg(544344918, ios == 0,
'unable to open file ' &
208 character(len=*),
intent(in) :: filename
210 integer,
intent(out) :: unit
215 open(unit=unit, file=filename, status=
'replace', action=
'write', &
217 call
assert_msg(609624199, ios == 0,
'unable to open file ' &
228 integer,
intent(in) :: unit
241 real(kind=dp),
intent(in) :: v
243 vol2rad = (v / (4d0 / 3d0 * const%pi))**(1d0/3d0)
253 real(kind=dp),
intent(in) :: v
255 vol2diam = 2d0 * (v / (4d0 / 3d0 * const%pi))**(1d0/3d0)
265 real(kind=dp),
intent(in) :: r
277 real(kind=dp),
intent(in) :: r
279 rad2vol = 4d0 / 3d0 * const%pi * r**3d0
289 real(kind=dp),
intent(in) :: d
301 real(kind=dp),
intent(in) :: d
303 diam2vol = 4d0 / 3d0 * const%pi * (d / 2d0)**3d0
314 real(kind=dp),
intent(in) :: d1
316 real(kind=dp),
intent(in) :: d2
319 real(kind=dp),
parameter :: eps = 1d-8
325 if (abs(d1 - d2) / (abs(d1) + abs(d2)) .lt. eps)
then
341 real(kind=dp),
intent(in) :: d1
343 real(kind=dp),
intent(in) :: d2
345 real(kind=dp),
intent(in) :: abs_tol
348 real(kind=dp),
parameter :: eps = 1d-8
354 if ((abs(d1 - d2) .lt. abs_tol) .or. &
355 (abs(d1 - d2) / (abs(d1) + abs(d2)) .lt. eps))
then
369 second_name, second_time)
372 character(len=*),
intent(in) :: first_name
374 real(kind=dp),
intent(in) :: first_time
376 character(len=*),
intent(in) :: second_name
378 real(kind=dp),
intent(in) :: second_time
380 real(kind=dp) :: ratio
382 ratio = first_time / second_time
383 if (abs(ratio - aint(ratio)) > 1d-6)
then
384 call
warn_msg(952299377, trim(first_name) &
385 //
" is not an integer multiple of " // trim(second_name))
403 real(kind=dp),
intent(in) :: time
405 real(kind=dp),
intent(in) :: timestep
407 real(kind=dp),
intent(in) :: interval
409 real(kind=dp),
intent(inout) :: last_time
411 logical,
intent(out) :: do_event
414 real(kind=dp),
parameter :: tolerance = 1d-6
416 real(kind=dp) closest_interval_time
419 if (time .eq. 0d0)
then
423 if ((time - last_time) .lt. interval * (1d0 - tolerance))
then
427 if ((time - last_time) .ge. interval)
then
432 closest_interval_time = anint(time / interval) * interval
433 if (abs(time - closest_interval_time) &
434 .lt. abs(time + timestep - closest_interval_time)) &
456 real(kind=dp),
intent(in) :: min_x
458 real(kind=dp),
intent(in) :: max_x
460 real(kind=dp),
intent(out) :: x(:)
467 a =
real(i - 1, kind=dp) /
real(n - 1, kind=dp)
468 x(i) = (1d0 - a) * min_x + a * max_x
484 real(kind=dp),
intent(in) :: min_x
486 real(kind=dp),
intent(in) :: max_x
488 real(kind=dp),
intent(out) :: x(:)
491 real(kind=dp) :: log_x(size(x))
495 call
assert(548290438, min_x > 0d0)
496 call
assert(805259035, max_x > 0d0)
497 call
linspace(log(min_x), log(max_x), log_x)
520 real(kind=dp),
intent(in) :: min_x
522 real(kind=dp),
intent(in) :: max_x
524 integer,
intent(in) :: n
526 real(kind=dp),
intent(in) :: x
533 *
real(n - 1, kind=dp)) + 1
554 real(kind=dp),
intent(in) :: min_x
556 real(kind=dp),
intent(in) :: max_x
558 integer,
intent(in) :: n
560 real(kind=dp),
intent(in) :: x
580 integer,
intent(in) :: n
582 real(kind=dp),
intent(in) :: x_vals(n)
584 real(kind=dp),
intent(in) :: x
593 do while (x >= x_vals(p))
614 real(kind=dp),
intent(in) :: x_vals(:)
616 real(kind=dp),
intent(in) :: y_vals(size(x_vals))
618 real(kind=dp),
intent(in) :: x
621 real(kind=dp) :: y, alpha
630 alpha = (x - x_vals(p)) / (x_vals(p+1) - x_vals(p))
631 y = (1d0 - alpha) * y_vals(p) + alpha * y_vals(p+1)
647 real(kind=dp),
intent(in) :: x_1
649 real(kind=dp),
intent(in) :: x_n
651 integer,
intent(in) :: n
653 integer,
intent(in) :: i
655 real(kind=dp) :: alpha
657 alpha =
real(i - 1, kind=dp) /
real(n - 1, kind=dp)
668 character(len=*),
intent(in) :: string
672 character(len=len(string)+300) :: error_msg
674 call
assert_msg(447772570, len_trim(string) <= 20, &
675 'error converting "' // trim(string) &
676 //
'" to integer: string too long')
677 read(string,
'(i20)', iostat=ios) val
679 'error converting "' // trim(string) &
691 character(len=*),
intent(in) :: string
695 character(len=len(string)+300) :: error_msg
697 call
assert_msg(733728030, len_trim(string) <= 30, &
698 'error converting "' // trim(string) //
'" to real: string too long')
699 read(string,
'(f30.0)', iostat=ios) val
701 'error converting "' // trim(string) &
713 character(len=*),
intent(in) :: string
717 character(len=len(string)+300) :: error_msg
720 if ((trim(string) ==
'yes') &
721 .or. (trim(string) ==
'y') &
722 .or. (trim(string) ==
'true') &
723 .or. (trim(string) ==
't') &
724 .or. (trim(string) ==
'1'))
then
726 elseif ((trim(string) ==
'no') &
727 .or. (trim(string) ==
'n') &
728 .or. (trim(string) ==
'false') &
729 .or. (trim(string) ==
'f') &
730 .or. (trim(string) ==
'0'))
then
733 call
die_msg(985010153,
'error converting "' // trim(string) &
746 integer,
intent(in) :: val
748 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
751 write(ret_val,
'(i30)') val
762 real(kind=dp),
intent(in) :: val
764 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
767 write(ret_val,
'(g30.20)') val
778 logical,
intent(in) :: val
780 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
798 complex(kind=dc),
intent(in) :: val
800 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
812 character(len=PMC_UTIL_CONVERT_STRING_LEN) &
816 integer,
intent(in) :: val
818 integer,
intent(in) :: max_len
820 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
823 if (len_trim(ret_val) > max_len)
then
833 character(len=PMC_UTIL_CONVERT_STRING_LEN) &
837 real(kind=dp),
intent(in) :: val
839 integer,
intent(in) :: max_len
841 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val, exp_str, frac_str
842 integer :: exp_val, exp_len, frac_len, use_frac_len, min_frac_len, i
843 real(kind=dp) :: frac_val
847 if (max_len >= 3)
then
858 exp_val = floor(log10(abs(val)))
859 frac_val = val / 10d0**exp_val
863 exp_len = len_trim(exp_str)
864 frac_len = len_trim(frac_str)
865 use_frac_len = max_len - 1 - exp_len
866 if (use_frac_len > frac_len)
then
867 use_frac_len = frac_len
874 if (use_frac_len < min_frac_len)
then
879 ret_val = frac_str(1:use_frac_len) //
"e" // trim(exp_str)
888 character(len=PMC_UTIL_CONVERT_STRING_LEN) &
892 real(kind=dp),
intent(in) :: time
894 integer,
intent(in) :: max_len
896 integer,
dimension(4),
parameter :: scale = (/ 1, 60, 60, 24 /)
897 character,
dimension(4),
parameter :: unit = (/
"s",
"m",
"h",
"d" /)
899 character(len=PMC_UTIL_CONVERT_STRING_LEN) :: ret_val
902 real(kind=dp) :: scaled_time
907 scaled_time = scaled_time /
real(scale(i), kind=dp)
909 if (len_trim(ret_val) <= max_len)
then
914 if (.not. len_ok)
then
933 integer,
intent(in) :: n
935 real(kind=dp),
intent(in) :: vec_cts(n)
937 integer,
intent(in) :: n_samp
939 integer,
intent(out) :: vec_disc(n)
942 real(kind=dp) :: vec_tot
944 vec_tot = sum(vec_cts)
947 vec_disc = nint(vec_cts / vec_tot *
real(n_samp, kind=dp))
950 do while (sum(vec_disc) < n_samp)
951 k = minloc(abs(
real(vec_disc + 1, kind=dp) - vec_cts) &
952 - abs(
real(vec_disc, kind=dp) - vec_cts))
953 vec_disc(k) = vec_disc(k) + 1
957 do while (sum(vec_disc) > n_samp)
958 k = minloc(abs(
real(vec_disc - 1, kind=dp) - vec_cts) &
959 - abs(
real(vec_disc, kind=dp) - vec_cts))
960 vec_disc(k) = vec_disc(k) - 1
963 call
assert_msg(323412496, sum(vec_disc) == n_samp, &
964 'generated incorrect number of samples')
974 integer,
intent(in) :: int_vec(:)
976 integer,
intent(out) :: int_avg
978 int_avg = sum(int_vec) /
size(int_vec)
988 real(kind=dp),
intent(in) :: real_vec(:)
990 real(kind=dp),
intent(out) :: real_avg
992 real_avg = sum(real_vec) /
real(size(real_vec), kind=dp)
1003 real(kind=dp),
intent(inout),
allocatable :: x(:)
1005 integer,
intent(in) :: n
1007 logical,
intent(in),
optional :: only_grow
1010 real(kind=dp),
allocatable :: tmp_x(:)
1012 if (
allocated(x))
then
1014 if (present(only_grow))
then
1015 new_n = max(new_n,
size(x))
1017 if (
size(x) /= new_n)
then
1018 allocate(tmp_x(new_n))
1020 tmp_x(1:min(new_n,
size(x))) = x(1:min(new_n,
size(x)))
1021 call move_alloc(tmp_x, x)
1037 real(kind=dp),
intent(inout),
allocatable :: x(:, :)
1039 integer,
intent(in) :: n1
1041 integer,
intent(in) :: n2
1043 logical,
intent(in),
optional :: only_grow
1045 integer :: new_n1, new_n2, n1_min, n2_min
1046 real(kind=dp),
allocatable :: tmp_x(:, :)
1048 if (
allocated(x))
then
1051 if (present(only_grow))
then
1052 new_n1 = max(new_n1,
size(x, 1))
1053 new_n2 = max(new_n2,
size(x, 2))
1055 if ((
size(x, 1) /= new_n1) .or. (
size(x, 2) /= new_n2))
then
1056 allocate(tmp_x(new_n1, new_n2))
1057 n1_min = min(new_n1,
size(x, 1))
1058 n2_min = min(new_n2,
size(x, 2))
1060 tmp_x(1:n1_min, 1:n2_min) = x(1:n1_min, 1:n2_min)
1061 call move_alloc(tmp_x, x)
1077 integer,
intent(inout),
allocatable :: x(:)
1079 integer,
intent(in) :: n
1081 logical,
intent(in),
optional :: only_grow
1084 integer,
allocatable :: tmp_x(:)
1086 if (
allocated(x))
then
1088 if (present(only_grow))
then
1089 new_n = max(new_n,
size(x))
1091 if (
size(x) /= new_n)
then
1092 allocate(tmp_x(new_n))
1094 tmp_x(1:min(new_n,
size(x))) = x(1:min(new_n,
size(x)))
1095 call move_alloc(tmp_x, x)
1111 integer,
intent(inout),
allocatable :: x(:, :)
1113 integer,
intent(in) :: n1
1115 integer,
intent(in) :: n2
1117 logical,
intent(in),
optional :: only_grow
1119 integer :: new_n1, new_n2, n1_min, n2_min
1120 integer,
allocatable :: tmp_x(:, :)
1122 if (
allocated(x))
then
1125 if (present(only_grow))
then
1126 new_n1 = max(new_n1,
size(x, 1))
1127 new_n2 = max(new_n2,
size(x, 2))
1129 if ((
size(x, 1) /= new_n1) .or. (
size(x, 2) /= new_n2))
then
1130 allocate(tmp_x(new_n1, new_n2))
1131 n1_min = min(new_n1,
size(x, 1))
1132 n2_min = min(new_n2,
size(x, 2))
1134 tmp_x(1:n1_min, 1:n2_min) = x(1:n1_min, 1:n2_min)
1135 call move_alloc(tmp_x, x)
1151 character(len=*),
intent(inout),
allocatable :: x(:)
1153 integer,
intent(in) :: n
1155 if (
allocated(x))
then
1156 if (
size(x) /= n)
then
1176 character(len=*),
intent(in) :: filename
1178 character(len=*),
intent(out) :: basename
1181 logical :: found_period
1184 i = len_trim(basename)
1185 found_period = .false.
1186 do while ((i > 0) .and. (.not. found_period))
1191 if (basename(i:i) ==
".")
then
1192 found_period = .true.
1209 character(len=*),
intent(out) :: date_time
1211 character(len=10) :: date, time, zone
1213 call
assert_msg(893219839, len(date_time) >= 29, &
1214 "date_time string must have length at least 29")
1215 call date_and_time(date, time, zone)
1217 write(date_time,
'(14a)') date(1:4),
"-", date(5:6),
"-", &
1218 date(7:8),
"T", time(1:2),
":", time(3:4),
":", &
1219 time(5:10), zone(1:3),
":", zone(4:5)
1229 real(kind=dp),
intent(in) :: deg
1231 deg2rad = deg / 180d0 * const%pi
1241 real(kind=dp),
intent(in) :: rad
1243 rad2deg = rad / const%pi * 180d0
1257 integer,
intent(inout) :: data(:)
1259 integer,
intent(out) :: perm(size(data))
1261 #ifdef PMC_USE_C_SORT
1262 integer(kind=c_int) :: n_c
1263 integer(kind=c_int),
target :: data_c(size(data))
1264 integer(kind=c_int),
target :: perm_c(size(data))
1265 type(c_ptr
) :: data_ptr, perm_ptr
1267 #ifndef DOXYGEN_SKIP_DOC
1271 integer(kind=c_int), value :: n_c
1272 type(c_ptr
), value :: data_ptr, perm_ptr
1277 data_c = int(
data, kind=c_int)
1279 n_c = int(
size(data), kind=c_int)
1280 data_ptr = c_loc(data_c)
1281 perm_ptr = c_loc(perm_c)
1295 character(len=*),
intent(in) :: filename
1297 real(kind=dp),
intent(inout),
allocatable :: data(:,:)
1299 integer :: unit, row, col
1301 character(len=1000) :: word
1310 do while (.not. eof)
1312 if (len_trim(word) > 0)
then
1314 if (col >
size(
data, 2))
then
1318 if (col >
size(
data, 2))
then
1319 call
assert_msg(516120334, col <=
size(
data, 2), &
1320 trim(filename) //
": line " &
1324 if (row >
size(
data, 1))
then
1330 if (eol .or. eof)
then
1335 (col == 1) .or. (col ==
size(
data, 2) + 1), &
1336 trim(filename) //
": line " &
1359 real(kind=dp),
allocatable,
intent(inout) :: data(:,:)
1361 integer,
intent(in) :: rows
1363 integer,
intent(in) :: cols
1365 real(kind=dp) :: tmp_data(rows, cols)
1366 integer :: data_rows, data_cols
1368 data_rows = min(rows,
size(
data, 1))
1369 data_cols = min(cols,
size(
data, 2))
1370 tmp_data(1:data_rows, 1:data_cols) =
data(1:data_rows, 1:data_cols)
1372 allocate(
data(rows, cols))
1373 data(1:data_rows, 1:data_cols) = tmp_data(1:data_rows, 1:data_cols)
1407 integer,
intent(in) :: unit
1409 character,
intent(out) :: char
1411 logical,
intent(out) :: eol
1413 logical,
intent(out) :: eof
1415 integer :: ios, n_read
1416 character(len=1) :: read_char
1422 read(unit=unit, fmt=
'(a)', advance=
'no', end=100, eor=110, &
1423 iostat=ios) read_char
1425 write(0,*)
'ERROR: reading file: IOSTAT = ', ios
1451 integer,
intent(in) :: unit
1453 character(len=*),
intent(out) :: word
1455 logical,
intent(out) :: eol
1457 logical,
intent(out) :: eof
1466 do while (((ichar(char) == 9) .or. (ichar(char) == 32)) &
1467 .and. (.not. eol) .and. (.not. eof))
1470 if (eol .or. eof)
return
1476 do while ((ichar(char) /= 9) .and. (ichar(char) /= 32) &
1477 .and. (.not. eol) .and. (.not. eof) .and. (i < len(word)))
1480 if (i < len(word))
then
1496 character(len=*),
intent(in) :: string
1498 character(len=*),
intent(in) :: start_string
1500 if (len(string) < len(start_string))
then
1504 if (string(1:len(start_string)) == start_string)
then
1518 real(kind=dp),
intent(in) :: source(:)
1520 real(kind=dp),
intent(inout),
pointer :: dest(:)
1522 if (
size(dest) /=
size(source))
then
1524 allocate(dest(
size(source)))
1536 real(kind=dp),
intent(in) :: source(:, :)
1538 real(kind=dp),
intent(inout),
pointer :: dest(:, :)
1540 if ((
size(dest, 1) /=
size(source, 1)) &
1541 .or. (
size(dest, 2) /=
size(source, 2)))
then
1543 allocate(dest(
size(source, 1),
size(source, 2)))
1555 real(kind=dp),
intent(in) :: x1
1557 real(kind=dp),
intent(in) :: x2
1569 real(kind=dp),
intent(in) :: p
1586 real(kind=dp),
intent(in) :: p(:)
real(kind=dp) elemental function vol2diam(v)
Convert volume (m^3) to diameter (m).
character(len=pmc_util_convert_string_len) function integer_to_string_max_len(val, max_len)
Convert an integer to a string format of maximum length.
subroutine read_char_raw(unit, char, eol, eof)
Read a single character from a file, signaling if we have hit EOL or EOF. If EOL or EOF are true then...
subroutine die_msg(code, error_msg)
Error immediately.
real(kind=dp) function rad2deg(rad)
Convert radians to degrees.
subroutine copy_real_2d(source, dest)
Copy a 2D array of reals, reallocating if necessary.
subroutine ensure_string_array_size(x, n)
Allocate or reallocate the given array to ensure it is of the given size, without preserving data...
integer function linspace_find(min_x, max_x, n, x)
Find the position of a real number in a 1D linear array.
real(kind=dp) elemental function diam2vol(d)
Convert diameter (m) to volume (m^3).
subroutine iso8601_date_and_time(date_time)
Current date and time in ISO 8601 format.
real(kind=dp) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
real(kind=dp) function entropy(p)
Compute the entropy of a probability mass function (non necessarily normalized).
int integer_sort_c(int n, int *data, int *perm)
Sort the given data array and return the permutation.
character(len=pmc_util_convert_string_len) function logical_to_string(val)
Convert a logical to a string format.
logical function starts_with(string, start_string)
Checks whether a string starts with a given other string.
subroutine close_file(unit)
Close a file and de-assign the unit.
integer function find_1d(n, x_vals, x)
Find the position of a real number in an arbitrary 1D array.
real(kind=dp) elemental function rad2diam(r)
Convert radius (m) to diameter (m).
subroutine ensure_integer_array_2d_size(x, n1, n2, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine logspace(min_x, max_x, x)
Makes a logarithmically spaced array of length n from min to max.
logical function almost_equal(d1, d2)
Tests whether two real numbers are almost equal using only a relative tolerance.
real(kind=dp) function interp_linear_disc(x_1, x_n, n, i)
Linear interpolation over discrete indices.
real(kind=dp) elemental function diam2rad(d)
Convert diameter (m) to radius (m).
integer function string_to_integer(string)
Convert a string to an integer.
subroutine check_time_multiple(first_name, first_time, second_name, second_time)
Check that the first time interval is close to an integer multiple of the second, and warn if it is n...
subroutine average_integer(int_vec, int_avg)
Computes the average of an array of integer numbers.
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
subroutine copy_real_1d(source, dest)
Copy a 1D array of reals, reallocating if necessary.
Common utility subroutines.
subroutine reallocate_real_array2d(data, rows, cols)
Reallocate a 2D real array to the given size, preserving the contents.
logical function string_to_logical(string)
Convert a string to a logical.
subroutine average_real(real_vec, real_avg)
Computes the average of an array of real numbers.
real(kind=dp) function deg2rad(deg)
Convert degrees to radians.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
character(len=pmc_util_convert_string_len) function complex_to_string(val)
Convert a complex to a string format.
logical function almost_equal_abs(d1, d2, abs_tol)
Tests whether two real numbers are almost equal using an absolute and relative tolerance.
subroutine die(code)
Error immediately.
integer function get_unit()
Returns an available unit number. This should be freed by free_unit().
subroutine ensure_real_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine warn_msg(code, warning_msg, already_warned)
Prints a warning message.
subroutine loadtxt(filename, data)
Load a real array from a text file.
subroutine check_event(time, timestep, interval, last_time, do_event)
Computes whether an event is scheduled to take place.
subroutine ensure_real_array_2d_size(x, n1, n2, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
elemental real(kind=dp) function harmonic_mean(x1, x2)
Compute the harmonic mean of two numbers.
subroutine integer_sort(data, perm)
Sort the given data array and return the permutation defining the sort.
real(kind=dp) function interp_1d(x_vals, y_vals, x)
1D linear interpolation.
real(kind=dp) elemental function vol2rad(v)
Convert volume (m^3) to radius (m).
character(len=pmc_util_convert_string_len) function real_to_string(val)
Convert a real to a string format.
integer function logspace_find(min_x, max_x, n, x)
Find the position of a real number in a 1D logarithmic array.
subroutine linspace(min_x, max_x, x)
Makes a linearly spaced array from min to max.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine free_unit(unit)
Frees a unit number returned by get_unit().
real(kind=dp) function string_to_real(string)
Convert a string to a real.
subroutine read_word_raw(unit, word, eol, eof)
Read a white-space delimited word from a file, signaling if we have EOL or EOF. If EOL or EOF are tru...
character(len=pmc_util_convert_string_len) function time_to_string_max_len(time, max_len)
Convert a time to a string format of maximum length.
subroutine vec_cts_to_disc(n, vec_cts, n_samp, vec_disc)
Convert a real-valued vector into an integer-valued vector.
character(len=pmc_util_convert_string_len) function real_to_string_max_len(val, max_len)
Convert a real to a string format of maximum length.
subroutine open_file_read(filename, unit)
Open a file for reading with an automatically assigned unit and test that it succeeds. The file should be closed with close_file().
elemental real(kind=dp) function nplogp(p)
Compute for computing entropy.
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine open_file_write(filename, unit)
Open a file for writing with an automatically assigned unit and test that it succeeds. The file should be closed with close_file().
subroutine get_basename(filename, basename)
Strip the extension to find the basename.