PLUME-MoM-TSM  1.0
VolcanicPlumeModel
complexify.f90
Go to the documentation of this file.
1 !******************************************************************************
2 ! Written for 'complexify.py 1.3'
3 ! J.R.R.A.Martins 1999
4 ! 21-Apr-00 Fixed tan, sinh, cosh
5 ! sign now returns complex
6 ! added log10 and nint
7 ! changed ==, /= and >= -- see comments below
8 ! 20-May-00 added cosd, sind, and epsilon
9 ! 11-Jul-00 took away cosd, sind (they are reserved, but not
10 ! intrinsic functions in F90)
11 ! 21-Jul-00 converted all trig functions to the value/derivative
12 ! formulas -- not general complex number formulas
13 ! 15-Aug-00 Fixed bug in atan2 formula and added the rest of the
14 ! _ci and _ic cominations to the relational operators.
15 ! P. Sturdza
16 !
17 !******************************************************************************
18 !
19 ! Assume all code is compiled with double precision (-r8 compiler flag)
20 !
21 
22 !TODO:
23 ! more typ combinations: cc, cr, rc, ic ?
24 ! check all fcns
25 !
26 
27 module complexify
28 
29  USE variables, ONLY : wp
30 
31  implicit none
32 
33 ! ABS
34  interface abs
35  module procedure abs_c
36  end interface
37 
38 ! COSD
39 ! interface cosd
40 ! module procedure cosd_c
41 ! end interface
42 
43 ! ACOS
44  interface acos
45  module procedure acos_c
46  end interface
47 
48 ! SIND
49 ! interface sind
50 ! module procedure sind_c
51 ! end interface
52 
53 ! ASIN
54  interface asin
55  module procedure asin_c
56  end interface
57 
58 ! ATAN
59  interface atan
60  module procedure atan_c
61  end interface
62 
63 ! ATAN2
64  interface atan2
65  module procedure atan2_cc
66  end interface
67 
68 ! COSH
69  interface cosh
70  module procedure cosh_c
71  end interface
72 
73 ! MAX (limited to 2-4 complex args, 2 mixed args)
74  interface max
75  module procedure max_cc
76  module procedure max_cr
77  module procedure max_rc
78  module procedure max_ccc ! added because of DFLUX.f
79  module procedure max_cccc ! added because of DFLUX.f
80  end interface
81 
82 ! MIN (limited to 2-4 complex args, 2 mixed args)
83  interface min
84  module procedure min_cc
85  module procedure min_cr
86  module procedure min_rc
87  module procedure min_ccc
88  module procedure min_cccc
89  end interface
90 
91 ! SIGN
92  interface sign
93  module procedure sign_cc
94  module procedure sign_cr
95  module procedure sign_rc
96  end interface
97 
98 ! DIM
99  interface dim
100  module procedure dim_cc
101  module procedure dim_cr
102  module procedure dim_rc
103  end interface
104 
105 ! SINH
106  interface sinh
107  module procedure sinh_c
108  end interface
109 
110 ! TAN
111  interface tan
112  module procedure tan_c
113  end interface
114 
115 ! TANH
116  interface tanh
117  module procedure tanh_c
118  end interface
119 
120 ! LOG10
121  interface log10
122  module procedure log10_c
123  end interface
124 
125 ! NINT
126  interface nint
127  module procedure nint_c
128  end interface
129 
130 ! EPSILON
131  interface epsilon
132  module procedure epsilon_c
133  end interface
134 
135 ! <
136  interface operator (<)
137  module procedure lt_cc
138  module procedure lt_cr
139  module procedure lt_rc
140  module procedure lt_ci
141  module procedure lt_ic
142  end interface
143 
144 ! <=
145  interface operator (<=)
146  module procedure le_cc
147  module procedure le_cr
148  module procedure le_rc
149  module procedure le_ci
150  module procedure le_ic
151  end interface
152 
153 ! >
154  interface operator (>)
155  module procedure gt_cc
156  module procedure gt_cr
157  module procedure gt_rc
158  module procedure gt_ci
159  module procedure gt_ic
160  end interface
161 
162 !! MIPSpro Compilers: Version 7.30 won't take .ge. and .eq..
163 !! But pgf90 on Linux doesn't complain, go figure.
164 !! It looks like a strict interpretation of FORTRAN should
165 !! not allow overloading of .eq. and .ne. since they already
166 !! have a definition for type complex, so define new operators
167 !! called .ceq., .cne. and, for MIPS, .cge.
168 !!
169 !! comment out (and uncomment) the appropriate versions for
170 !! your compiler
171 !!
172 ! >=
173  interface operator (>=)
174  module procedure ge_cc
175  module procedure ge_cr
176  module procedure ge_rc
177  module procedure ge_ci
178  module procedure ge_ic
179  end interface
180 ! interface operator (.cge.)
181 ! module procedure ge_cc
182 ! module procedure ge_rr
183 ! module procedure ge_ii
184 ! module procedure ge_aa
185 ! module procedure ge_cr
186 ! module procedure ge_rc
187 ! module procedure ge_ci
188 ! module procedure ge_ic
189 ! module procedure ge_ir
190 ! module procedure ge_ri
191 ! end interface
192 
193 ! ==
194 ! interface operator (==)
195 ! module procedure eq_cc
196 ! module procedure eq_cr
197 ! module procedure eq_rc
198 ! module procedure eq_ci
199 ! module procedure eq_ic
200 ! end interface
201  interface operator (.ceq.)
202  module procedure eq_cc
203  module procedure eq_rr
204  module procedure eq_ii
205  module procedure eq_aa
206  module procedure eq_cr
207  module procedure eq_rc
208  module procedure eq_ci
209  module procedure eq_ic
210  module procedure eq_ir
211  module procedure eq_ri
212  end interface
213 
214 ! /=
215 ! interface operator (/=)
216 ! module procedure ne_cc
217 ! module procedure ne_cr
218 ! module procedure ne_rc
219 ! module procedure ne_ci
220 ! module procedure ne_ic
221 ! end interface
222  interface operator (.cne.)
223  module procedure ne_cc
224  module procedure ne_rr
225  module procedure ne_ii
226  module procedure ne_aa
227  module procedure ne_cr
228  module procedure ne_rc
229  module procedure ne_ci
230  module procedure ne_ic
231  module procedure ne_ir
232  module procedure ne_ri
233  end interface
234 
235 contains
236 
237 !******************************************************************************
238 !
239 ! Function definitions
240 !
241 !******************************************************************************
242 
243 ! ABS, intrinsic
244  complex*16 function abs_c(val)
245  complex*16, intent(in) :: val
246  abs_c = val
247  if (real(val) < 0) abs_c = cmplx(-real(val),-aimag(val))
248  return
249  end function abs_c
250 
251 ! COSD
252 ! complex*16 function cosd_c(z)
253 ! complex*16, intent(in) :: z
254 ! cosd_c = cos(z*3.14159265358979323846/180.)
255 ! end function cosd_c
256 
257 ! SIND
258 ! complex*16 function sind_c(z)
259 ! complex*16, intent(in) :: z
260 ! sind_c = sin(z*3.14159265358979323846/180.)
261 ! end function sind_c
262 
263 ! ACOS
264  complex*16 function acos_c(z)
265  complex*16, intent(in) :: z
266 ! acos_c = - cmplx(0., 1.)*log(z+sqrt(z**2-1.))
267 ! not general complex valued formula:
268  acos_c = cmplx(acos(real(z)),-aimag(z)/sqrt(1.-real(z)**2))
269  return
270  end function acos_c
271 
272 ! ASIN
273  complex*16 function asin_c(z)
274  complex*16, intent(in) :: z
275 ! asin_c = - cmplx(0., 1.)*log(cmplx(0.,1.)*z+sqrt(1.-z**2))
276 ! not general complex valued formula:
277  asin_c = cmplx(asin(real(z)),aimag(z)/sqrt(1.-real(z)**2))
278  return
279  end function asin_c
280 
281 ! ATAN
282  complex*16 function atan_c(z)
283  complex*16, intent(in) :: z
284 ! complex*16 z2
285 ! real(wp) pi2, xans, yans, r, r2, x, y
286 ! pi2 = 2.0*atan(1.0)
287 ! r = sqrt(real(z)**2+aimag(z)**2)
288 ! x = real(z)
289 ! y = aimag(z)
290 ! r2 = r*r
291 ! xans = 0.5*atan2 (2.0*x, 1.0-r2)
292 ! yans = 0.25*log((r2+2.0*y+1.0)/(r2-2.0*y+1.0))
293 ! atan_c = cmplx (xans, yans)
294 ! not general complex valued formula:
295  atan_c = cmplx(atan(real(z)),aimag(z)/(1.+real(z)**2))
296  return
297  end function atan_c
298 
299 ! ATAN2
300  complex*16 function atan2_cc(csn, ccs)
301  complex*16, intent(in) :: csn, ccs
302 ! real(wp) pi
303 ! pi = 4.0*atan(1.0)
304 ! if (sqrt(real(ccs)**2 + aimag(ccs)**2).eq.0.) then ! abs orig
305 ! if (sqrt(real(csn)**2+aimag(csn)**2).eq.0.) then
306 ! atan2_cc = cmplx(0.0)
307 ! else
308 ! atan2_cc = cmplx(sign(0.5*pi,real(csn)), 0.0)
309 ! end if
310 ! else
311 ! atan2_cc = atan(csn/ccs)
312 ! if (real(ccs).lt.0.) atan2_cc = atan2_cc + pi
313 ! if (real(atan2_cc).gt.pi) atan2_cc = atan2_cc - 2.0*pi
314 ! end if
315 ! not general complex valued formula:
316  real(wp) a,b,c,d
317  a=real(csn)
318  b=aimag(csn)
319  c=real(ccs)
320  d=aimag(ccs)
321  atan2_cc=cmplx(atan2(a,c),(c*b-a*d)/(a**2+c**2))
322  return
323  end function atan2_cc
324 
325 ! COSH
326  complex*16 function cosh_c(z)
327  complex*16, intent(in) :: z
328 ! complex*16 eplus, eminus
329 ! eplus = exp(z)
330 ! eminus = exp(z)
331 ! cosh_c = (eplus + eminus)/2.
332 ! not general complex valued formula:
333  cosh_c=cmplx(cosh(real(z)),aimag(z)*sinh(real(z)))
334  return
335  end function cosh_c
336 
337 ! SINH
338  complex*16 function sinh_c(z)
339  complex*16, intent(in) :: z
340 ! complex*16 eplus, eminus
341 ! eplus = exp(z)
342 ! eminus = exp(z)
343 ! sinh_c = (eplus - eminus)/2.
344 ! not general complex valued formula:
345  sinh_c=cmplx(sinh(real(z)),aimag(z)*cosh(real(z)))
346  return
347  end function sinh_c
348 
349 ! TAN
350  complex*16 function tan_c(z)
351  complex*16, intent(in) :: z
352 ! complex*16 eiplus, eiminus
353 ! eiplus = exp(cmplx(0.,1.)*z)
354 ! eiminus = exp(-cmplx(0.,1.)*z)
355 ! tan_c = cmplx(0.,1.)*(eiminus - eiplus)/(eiplus + eiminus)
356 ! not general complex valued formula:
357  tan_c=cmplx(tan(real(z)),aimag(z)/cos(real(z))**2)
358  return
359  end function tan_c
360 
361 ! TANH
362  complex*16 function tanh_c(a)
363  complex*16, intent(in) :: a
364 ! complex*16 eplus, eminus
365 ! if(real(a) > 50)then
366 ! tanh_c = 1.
367 ! else
368 ! eplus = exp(a)
369 ! eminus = exp(-a)
370 ! tanh_c = (eplus - eminus)/(eplus + eminus)
371 ! end if
372 ! not general complex valued formula:
373  tanh_c=cmplx(tanh(real(a)),aimag(a)/cosh(real(a))**2)
374  return
375  end function tanh_c
376 
377 ! MAX, intrinsic
378  complex*16 function max_cc(val1, val2)
379  complex*16, intent(in) :: val1, val2
380  if (real(val1) > real(val2)) then
381  max_cc = val1
382  else
383  max_cc = val2
384  endif
385  return
386  end function max_cc
387  complex*16 function max_cr(val1, val2)
388  complex*16, intent(in) :: val1
389  real(wp), intent(in) :: val2
390  if (real(val1) > val2) then
391  max_cr = val1
392  else
393  max_cr = cmplx(val2, 0.)
394  endif
395  return
396  end function max_cr
397  complex*16 function max_rc(val1, val2)
398  real(wp), intent(in) :: val1
399  complex*16, intent(in) :: val2
400  if (val1 > real(val2)) then
401  max_rc = cmplx(val1, 0.)
402  else
403  max_rc = val2
404  endif
405  return
406  end function max_rc
407  complex*16 function max_ccc(val1, val2, val3)
408  complex*16, intent(in) :: val1, val2, val3
409  if (real(val1) > real(val2)) then
410  max_ccc = val1
411  else
412  max_ccc = val2
413  endif
414  if (real(val3) > real(max_ccc)) then
415  max_ccc = val3
416  endif
417  return
418  end function max_ccc
419  function max_cccc(val1, val2, val3, val4)
420  complex*16, intent(in) :: val1, val2, val3, val4
421  complex*16 max_cccc
422  complex*16 max_cccc2
423  if (real(val1) > real(val2)) then
424  max_cccc = val1
425  else
426  max_cccc = val2
427  endif
428  if (real(val3) > real(val4)) then
429  max_cccc2 = val3
430  else
431  max_cccc2 = val4
432  endif
433  if (real(max_cccc2) > real(max_cccc)) then
434  max_cccc = max_cccc2
435  endif
436  return
437  end function max_cccc
438 
439 ! MIN, intrinsic
440  complex*16 function min_cc(val1, val2)
441  complex*16, intent(in) :: val1, val2
442  if (real(val1) < real(val2)) then
443  min_cc = val1
444  else
445  min_cc = val2
446  endif
447  return
448  end function min_cc
449  complex*16 function min_cr(val1, val2)
450  complex*16, intent(in) :: val1
451  real(wp), intent(in) :: val2
452  if (real(val1) < val2) then
453  min_cr = val1
454  else
455  min_cr = cmplx(val2, 0.)
456  endif
457  return
458  end function min_cr
459  complex*16 function min_rc(val1, val2)
460  real(wp), intent(in) :: val1
461  complex*16, intent(in) :: val2
462  if (val1 < real(val2)) then
463  min_rc = cmplx(val1, 0.)
464  else
465  min_rc = val2
466  endif
467  return
468  end function min_rc
469  complex*16 function min_ccc(val1, val2, val3)
470  complex*16, intent(in) :: val1, val2, val3
471  if (real(val1) < real(val2)) then
472  min_ccc = val1
473  else
474  min_ccc = val2
475  endif
476  if (real(val3) < real(min_ccc)) then
477  min_ccc = val3
478  endif
479  return
480  end function min_ccc
481  function min_cccc(val1, val2, val3, val4)
482  complex*16, intent(in) :: val1, val2, val3, val4
483  complex*16 min_cccc
484  complex*16 min_cccc2
485  if (real(val1) < real(val2)) then
486  min_cccc = val1
487  else
488  min_cccc = val2
489  endif
490  if (real(val3) < real(val4)) then
491  min_cccc2 = val3
492  else
493  min_cccc2 = val4
494  endif
495  if (real(min_cccc2) < real(min_cccc)) then
496  min_cccc = min_cccc2
497  endif
498  return
499  end function min_cccc
500 
501 
502 ! SIGN, intrinsic, assume that val1 is always a complex*16
503 ! in reality could be int
504  complex*16 function sign_cc(val1, val2)
505  complex*16, intent(in) :: val1, val2
506  real(wp) sign
507  if (real(val2) < 0.) then
508  sign = -1.
509  else
510  sign = 1.
511  endif
512  sign_cc = sign * val1
513  return
514  end function sign_cc
515  complex*16 function sign_cr(val1, val2)
516  complex*16, intent(in) :: val1
517  real(wp), intent(in) :: val2
518  real(wp) sign
519  if (real(val2) < 0.) then
520  sign = -1.
521  else
522  sign = 1.
523  endif
524  sign_cr = sign * val1
525  return
526  end function sign_cr
527  complex*16 function sign_rc(val1, val2)
528  real(wp), intent(in) :: val1
529  complex*16, intent(in) :: val2
530  real(wp) sign
531  if (real(val2) < 0.) then
532  sign = -1.
533  else
534  sign = 1.
535  endif
536  sign_rc = sign * val1
537  return
538  end function sign_rc
539 
540 ! DIM, intrinsic
541  complex*16 function dim_cc(val1, val2)
542  complex*16, intent(in) :: val1, val2
543  if (val1 > val2) then
544  dim_cc = val1 - val2
545  else
546  dim_cc = cmplx(0., 0.)
547  endif
548  return
549  end function dim_cc
550  complex*16 function dim_cr(val1, val2)
551  complex*16, intent(in) :: val1
552  real(wp), intent(in) :: val2
553  if (val1 > val2) then
554  dim_cr = val1 - cmplx(val2, 0.)
555  else
556  dim_cr = cmplx(0., 0.)
557  endif
558  return
559  end function dim_cr
560  complex*16 function dim_rc(val1, val2)
561  real(wp), intent(in) :: val1
562  complex*16, intent(in) :: val2
563  if (val1 > val2) then
564  dim_rc = cmplx(val1, 0.) - val2
565  else
566  dim_rc = cmplx(0., 0.)
567  endif
568  return
569  end function dim_rc
570 
571 ! LOG10
572  complex*16 function log10_c(z)
573  complex*16, intent(in) :: z
574  log10_c=log(z)/log((10.0,0.0))
575  end function log10_c
576 
577 ! NINT
578  integer function nint_c(z)
579  complex*16, intent(in) :: z
580  nint_c=nint(real(z))
581  end function nint_c
582 
583 ! EPSILON !! bad news ulness compiled with -r8
584  complex*16 function epsilon_c(z)
585  complex*16, intent(in) :: z
586  epsilon_c=epsilon(real(z))
587  end function epsilon_c
588 
589 ! <, .lt.
590  logical function lt_cc(lhs, rhs)
591  complex*16, intent(in) :: lhs, rhs
592  lt_cc = real(lhs) < real(rhs)
593  end function lt_cc
594  logical function lt_cr(lhs, rhs)
595  complex*16, intent(in) :: lhs
596  real(wp), intent(in) :: rhs
597  lt_cr = real(lhs) < rhs
598  end function lt_cr
599  logical function lt_rc(lhs, rhs)
600  real(wp), intent(in) :: lhs
601  complex*16, intent(in) :: rhs
602  lt_rc = lhs < real(rhs)
603  end function lt_rc
604  logical function lt_ci(lhs, rhs)
605  complex*16, intent(in) :: lhs
606  integer, intent(in) :: rhs
607  lt_ci = real(lhs) < rhs
608  end function lt_ci
609  logical function lt_ic(lhs, rhs)
610  integer, intent(in) :: lhs
611  complex*16, intent(in) :: rhs
612  lt_ic = lhs < real(rhs)
613  end function lt_ic
614 
615 ! <=, .le.
616  logical function le_cc(lhs, rhs)
617  complex*16, intent(in) :: lhs, rhs
618  le_cc = real(lhs) <= real(rhs)
619  end function le_cc
620  logical function le_cr(lhs, rhs)
621  complex*16, intent(in) :: lhs
622  real(wp), intent(in) :: rhs
623  le_cr = real(lhs) <= rhs
624  end function le_cr
625  logical function le_rc(lhs, rhs)
626  real(wp), intent(in) :: lhs
627  complex*16, intent(in) :: rhs
628  le_rc = lhs <= real(rhs)
629  end function le_rc
630  logical function le_ci(lhs, rhs)
631  complex*16, intent(in) :: lhs
632  integer, intent(in) :: rhs
633  le_ci = real(lhs) <= rhs
634  end function le_ci
635  logical function le_ic(lhs, rhs)
636  integer, intent(in) :: lhs
637  complex*16, intent(in) :: rhs
638  le_ic = lhs <= real(rhs)
639  end function le_ic
640 
641 ! >, .gt.
642  logical function gt_cc(lhs, rhs)
643  complex*16, intent(in) :: lhs, rhs
644  gt_cc = real(lhs) > real(rhs)
645  end function gt_cc
646  logical function gt_cr(lhs, rhs)
647  complex*16, intent(in) :: lhs
648  real(wp), intent(in) :: rhs
649  gt_cr = real(lhs) > rhs
650  end function gt_cr
651  logical function gt_rc(lhs, rhs)
652  real(wp), intent(in) :: lhs
653  complex*16, intent(in) :: rhs
654  gt_rc = lhs > real(rhs)
655  end function gt_rc
656  logical function gt_ci(lhs, rhs)
657  complex*16, intent(in) :: lhs
658  integer, intent(in) :: rhs
659  gt_ci = real(lhs) > rhs
660  end function gt_ci
661  logical function gt_ic(lhs, rhs)
662  integer, intent(in) :: lhs
663  complex*16, intent(in) :: rhs
664  gt_ic = lhs > real(rhs)
665  end function gt_ic
666 
667 !! here are the redefined ones:
668 ! >=, .ge.
669  logical function ge_cc(lhs, rhs)
670  complex*16, intent(in) :: lhs, rhs
671  ge_cc = real(lhs) >= real(rhs)
672  end function ge_cc
673  logical function ge_rr(lhs, rhs)
674  real(wp), intent(in) :: lhs, rhs
675  ge_rr = lhs >= rhs
676  end function ge_rr
677  logical function ge_ii(lhs, rhs)
678  integer, intent(in) :: lhs, rhs
679  ge_ii = lhs >= rhs
680  end function ge_ii
681  logical function ge_aa(lhs, rhs)
682  character(len=*), intent(in) :: lhs, rhs
683  ge_aa = lhs >= rhs
684  end function ge_aa
685  logical function ge_cr(lhs, rhs)
686  complex*16, intent(in) :: lhs
687  real(wp), intent(in) :: rhs
688  ge_cr = real(lhs) >= rhs
689  end function ge_cr
690  logical function ge_rc(lhs, rhs)
691  real(wp), intent(in) :: lhs
692  complex*16, intent(in) :: rhs
693  ge_rc = lhs >= real(rhs)
694  end function ge_rc
695  logical function ge_ci(lhs, rhs)
696  complex*16, intent(in) :: lhs
697  integer, intent(in) :: rhs
698  ge_ci = real(lhs) >= rhs
699  end function ge_ci
700  logical function ge_ic(lhs, rhs)
701  integer, intent(in) :: lhs
702  complex*16, intent(in) :: rhs
703  ge_ic = lhs >= real(rhs)
704  end function ge_ic
705  logical function ge_ir(lhs, rhs)
706  integer, intent(in) :: lhs
707  real(wp), intent(in) :: rhs
708  ge_ir = lhs >= rhs
709  end function ge_ir
710  logical function ge_ri(lhs, rhs)
711  real(wp), intent(in) :: lhs
712  integer, intent(in) :: rhs
713  ge_ri = lhs >= rhs
714  end function ge_ri
715 
716 ! ==, .eq.
717  logical function eq_cc(lhs, rhs)
718  complex*16, intent(in) :: lhs, rhs
719  eq_cc = real(lhs) == real(rhs)
720  end function eq_cc
721  logical function eq_rr(lhs, rhs)
722  real(wp), intent(in) :: lhs, rhs
723  eq_rr = lhs == rhs
724  end function eq_rr
725  logical function eq_ii(lhs, rhs)
726  integer, intent(in) :: lhs, rhs
727  eq_ii = lhs == rhs
728  end function eq_ii
729  logical function eq_aa(lhs, rhs)
730  character(len=*), intent(in) :: lhs, rhs
731  eq_aa = lhs == rhs
732  end function eq_aa
733  logical function eq_cr(lhs, rhs)
734  complex*16, intent(in) :: lhs
735  real(wp), intent(in) :: rhs
736  eq_cr = real(lhs) == rhs
737  end function eq_cr
738  logical function eq_rc(lhs, rhs)
739  real(wp), intent(in) :: lhs
740  complex*16, intent(in) :: rhs
741  eq_rc = lhs == real(rhs)
742  end function eq_rc
743  logical function eq_ci(lhs, rhs)
744  complex*16, intent(in) :: lhs
745  integer, intent(in) :: rhs
746  eq_ci = real(lhs) == rhs
747  end function eq_ci
748  logical function eq_ic(lhs, rhs)
749  integer, intent(in) :: lhs
750  complex*16, intent(in) :: rhs
751  eq_ic = lhs == real(rhs)
752  end function eq_ic
753  logical function eq_ir(lhs, rhs)
754  integer, intent(in) :: lhs
755  real(wp), intent(in) :: rhs
756  eq_ir = lhs == rhs
757  end function eq_ir
758  logical function eq_ri(lhs, rhs)
759  real(wp), intent(in) :: lhs
760  integer, intent(in) :: rhs
761  eq_ri = lhs == rhs
762  end function eq_ri
763 
764 ! /=, .ne.
765  logical function ne_cc(lhs, rhs)
766  complex*16, intent(in) :: lhs, rhs
767  ne_cc = real(lhs) /= real(rhs)
768  end function ne_cc
769  logical function ne_rr(lhs, rhs)
770  real(wp), intent(in) :: lhs, rhs
771  ne_rr = lhs /= rhs
772  end function ne_rr
773  logical function ne_ii(lhs, rhs)
774  integer, intent(in) :: lhs, rhs
775  ne_ii = lhs /= rhs
776  end function ne_ii
777  logical function ne_aa(lhs, rhs)
778  character(len=*), intent(in) :: lhs, rhs
779  ne_aa = lhs /= rhs
780  end function ne_aa
781  logical function ne_cr(lhs, rhs)
782  complex*16, intent(in) :: lhs
783  real(wp), intent(in) :: rhs
784  ne_cr = real(lhs) /= rhs
785  end function ne_cr
786  logical function ne_rc(lhs, rhs)
787  real(wp), intent(in) :: lhs
788  complex*16, intent(in) :: rhs
789  ne_rc = lhs /= real(rhs)
790  end function ne_rc
791  logical function ne_ci(lhs, rhs)
792  complex*16, intent(in) :: lhs
793  integer, intent(in) :: rhs
794  ne_ci = real(lhs) /= rhs
795  end function ne_ci
796  logical function ne_ic(lhs, rhs)
797  integer, intent(in) :: lhs
798  complex*16, intent(in) :: rhs
799  ne_ic = lhs /= real(rhs)
800  end function ne_ic
801  logical function ne_ir(lhs, rhs)
802  integer, intent(in) :: lhs
803  real(wp), intent(in) :: rhs
804  ne_ir = lhs /= rhs
805  end function ne_ir
806  logical function ne_ri(lhs, rhs)
807  real(wp), intent(in) :: lhs
808  integer, intent(in) :: rhs
809  ne_ri = lhs /= rhs
810  end function ne_ri
811 
812 end module complexify
complex *16 function max_ccc(val1, val2, val3)
Definition: complexify.f90:408
logical function ge_cc(lhs, rhs)
Definition: complexify.f90:670
logical function eq_ii(lhs, rhs)
Definition: complexify.f90:726
logical function le_ic(lhs, rhs)
Definition: complexify.f90:636
logical function eq_cc(lhs, rhs)
Definition: complexify.f90:718
complex *16 function max_cc(val1, val2)
Definition: complexify.f90:379
logical function lt_rc(lhs, rhs)
Definition: complexify.f90:600
complex *16 function sinh_c(z)
Definition: complexify.f90:339
logical function eq_ci(lhs, rhs)
Definition: complexify.f90:744
logical function ne_ii(lhs, rhs)
Definition: complexify.f90:774
logical function ne_aa(lhs, rhs)
Definition: complexify.f90:778
complex *16 function atan2_cc(csn, ccs)
Definition: complexify.f90:301
logical function lt_ic(lhs, rhs)
Definition: complexify.f90:610
logical function ge_ir(lhs, rhs)
Definition: complexify.f90:706
logical function gt_ic(lhs, rhs)
Definition: complexify.f90:662
complex *16 function min_cr(val1, val2)
Definition: complexify.f90:450
logical function ge_aa(lhs, rhs)
Definition: complexify.f90:682
complex *16 function min_rc(val1, val2)
Definition: complexify.f90:460
logical function gt_cc(lhs, rhs)
Definition: complexify.f90:643
complex *16 function dim_cr(val1, val2)
Definition: complexify.f90:551
logical function eq_rc(lhs, rhs)
Definition: complexify.f90:739
integer function nint_c(z)
Definition: complexify.f90:579
logical function ne_ic(lhs, rhs)
Definition: complexify.f90:797
complex *16 function max_cr(val1, val2)
Definition: complexify.f90:388
complex *16 function sign_rc(val1, val2)
Definition: complexify.f90:528
complex *16 function min_cccc(val1, val2, val3, val4)
Definition: complexify.f90:482
logical function eq_ic(lhs, rhs)
Definition: complexify.f90:749
logical function le_ci(lhs, rhs)
Definition: complexify.f90:631
logical function gt_cr(lhs, rhs)
Definition: complexify.f90:647
complex *16 function tanh_c(a)
Definition: complexify.f90:363
logical function eq_ri(lhs, rhs)
Definition: complexify.f90:759
complex *16 function log10_c(z)
Definition: complexify.f90:573
complex *16 function dim_cc(val1, val2)
Definition: complexify.f90:542
logical function ge_cr(lhs, rhs)
Definition: complexify.f90:686
logical function ge_ri(lhs, rhs)
Definition: complexify.f90:711
logical function ne_cc(lhs, rhs)
Definition: complexify.f90:766
logical function ne_ci(lhs, rhs)
Definition: complexify.f90:792
complex *16 function sign_cc(val1, val2)
Definition: complexify.f90:505
logical function eq_rr(lhs, rhs)
Definition: complexify.f90:722
complex *16 function acos_c(z)
Definition: complexify.f90:265
logical function gt_rc(lhs, rhs)
Definition: complexify.f90:652
complex *16 function sign_cr(val1, val2)
Definition: complexify.f90:516
logical function ge_ic(lhs, rhs)
Definition: complexify.f90:701
logical function ne_rc(lhs, rhs)
Definition: complexify.f90:787
complex *16 function min_cc(val1, val2)
Definition: complexify.f90:441
logical function ne_cr(lhs, rhs)
Definition: complexify.f90:782
logical function ge_rr(lhs, rhs)
Definition: complexify.f90:674
logical function ge_ci(lhs, rhs)
Definition: complexify.f90:696
logical function ne_ri(lhs, rhs)
Definition: complexify.f90:807
logical function le_cr(lhs, rhs)
Definition: complexify.f90:621
logical function ne_rr(lhs, rhs)
Definition: complexify.f90:770
complex *16 function asin_c(z)
Definition: complexify.f90:274
logical function ne_ir(lhs, rhs)
Definition: complexify.f90:802
logical function ge_ii(lhs, rhs)
Definition: complexify.f90:678
logical function lt_cc(lhs, rhs)
Definition: complexify.f90:591
logical function le_rc(lhs, rhs)
Definition: complexify.f90:626
integer, parameter wp
working precision
Definition: variables.f90:21
complex *16 function max_rc(val1, val2)
Definition: complexify.f90:398
logical function gt_ci(lhs, rhs)
Definition: complexify.f90:657
complex *16 function min_ccc(val1, val2, val3)
Definition: complexify.f90:470
logical function ge_rc(lhs, rhs)
Definition: complexify.f90:691
complex *16 function abs_c(val)
Definition: complexify.f90:245
logical function le_cc(lhs, rhs)
Definition: complexify.f90:617
complex *16 function cosh_c(z)
Definition: complexify.f90:327
complex *16 function atan_c(z)
Definition: complexify.f90:283
logical function eq_ir(lhs, rhs)
Definition: complexify.f90:754
complex *16 function dim_rc(val1, val2)
Definition: complexify.f90:561
complex *16 function epsilon_c(z)
Definition: complexify.f90:585
complex *16 function tan_c(z)
Definition: complexify.f90:351
complex *16 function max_cccc(val1, val2, val3, val4)
Definition: complexify.f90:420
Global variables.
Definition: variables.f90:10
logical function eq_aa(lhs, rhs)
Definition: complexify.f90:730
logical function lt_ci(lhs, rhs)
Definition: complexify.f90:605
logical function lt_cr(lhs, rhs)
Definition: complexify.f90:595
logical function eq_cr(lhs, rhs)
Definition: complexify.f90:734