chfie.f 3.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108
  1. *DECK CHFIE
  2. REAL FUNCTION CHFIE (X1, X2, F1, F2, D1, D2, A, B)
  3. C***BEGIN PROLOGUE CHFIE
  4. C***SUBSIDIARY
  5. C***PURPOSE Evaluates integral of a single cubic for PCHIA
  6. C***LIBRARY SLATEC (PCHIP)
  7. C***TYPE SINGLE PRECISION (CHFIE-S, DCHFIE-D)
  8. C***AUTHOR Fritsch, F. N., (LLNL)
  9. C***DESCRIPTION
  10. C
  11. C CHFIE: Cubic Hermite Function Integral Evaluator.
  12. C
  13. C Called by PCHIA to evaluate the integral of a single cubic (in
  14. C Hermite form) over an arbitrary interval (A,B).
  15. C
  16. C ----------------------------------------------------------------------
  17. C
  18. C Calling sequence:
  19. C
  20. C REAL X1, X2, F1, F2, D1, D2, A, B
  21. C REAL VALUE, CHFIE
  22. C
  23. C VALUE = CHFIE (X1, X2, F1, F2, D1, D2, A, B)
  24. C
  25. C Parameters:
  26. C
  27. C VALUE -- (output) value of the requested integral.
  28. C
  29. C X1,X2 -- (input) endpoints if interval of definition of cubic.
  30. C
  31. C F1,F2 -- (input) function values at the ends of the interval.
  32. C
  33. C D1,D2 -- (input) derivative values at the ends of the interval.
  34. C
  35. C A,B -- (input) endpoints of interval of integration.
  36. C
  37. C***SEE ALSO PCHIA
  38. C***ROUTINES CALLED (NONE)
  39. C***REVISION HISTORY (YYMMDD)
  40. C 820730 DATE WRITTEN
  41. C 820805 Converted to SLATEC library version.
  42. C 870813 Minor cosmetic changes.
  43. C 890411 1. Added SAVE statements (Vers. 3.2).
  44. C 2. Added SIX to REAL declaration.
  45. C 890411 REVISION DATE from Version 3.2
  46. C 891214 Prologue converted to Version 4.0 format. (BAB)
  47. C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
  48. C 900328 Added TYPE section. (WRB)
  49. C 910408 Updated AUTHOR section in prologue. (WRB)
  50. C 930503 Corrected to set VALUE=0 when IERR.ne.0. (FNF)
  51. C 930504 Eliminated IERR and changed name from CHFIV to CHFIE. (FNF)
  52. C***END PROLOGUE CHFIE
  53. C
  54. C Programming notes:
  55. C 1. There is no error return from this routine because zero is
  56. C indeed the mathematically correct answer when X1.EQ.X2 .
  57. C**End
  58. C
  59. C DECLARE ARGUMENTS.
  60. C
  61. REAL X1, X2, F1, F2, D1, D2, A, B
  62. C
  63. C DECLARE LOCAL VARIABLES.
  64. C
  65. REAL DTERM, FOUR, FTERM, H, HALF, PHIA1, PHIA2, PHIB1, PHIB2,
  66. * PSIA1, PSIA2, PSIB1, PSIB2, SIX, TA1, TA2, TB1, TB2, THREE,
  67. * TWO, UA1, UA2, UB1, UB2
  68. SAVE HALF, TWO, THREE, FOUR, SIX
  69. C
  70. C INITIALIZE.
  71. C
  72. DATA HALF /0.5/, TWO /2./, THREE /3./, FOUR /4./, SIX /6./
  73. C
  74. C VALIDITY CHECK INPUT.
  75. C
  76. C***FIRST EXECUTABLE STATEMENT CHFIE
  77. IF (X1 .EQ. X2) THEN
  78. CHFIE = 0
  79. ELSE
  80. H = X2 - X1
  81. TA1 = (A - X1) / H
  82. TA2 = (X2 - A) / H
  83. TB1 = (B - X1) / H
  84. TB2 = (X2 - B) / H
  85. C
  86. UA1 = TA1**3
  87. PHIA1 = UA1 * (TWO - TA1)
  88. PSIA1 = UA1 * (THREE*TA1 - FOUR)
  89. UA2 = TA2**3
  90. PHIA2 = UA2 * (TWO - TA2)
  91. PSIA2 = -UA2 * (THREE*TA2 - FOUR)
  92. C
  93. UB1 = TB1**3
  94. PHIB1 = UB1 * (TWO - TB1)
  95. PSIB1 = UB1 * (THREE*TB1 - FOUR)
  96. UB2 = TB2**3
  97. PHIB2 = UB2 * (TWO - TB2)
  98. PSIB2 = -UB2 * (THREE*TB2 - FOUR)
  99. C
  100. FTERM = F1*(PHIA2 - PHIB2) + F2*(PHIB1 - PHIA1)
  101. DTERM = ( D1*(PSIA2 - PSIB2) + D2*(PSIB1 - PSIA1) )*(H/SIX)
  102. C
  103. CHFIE = (HALF*H) * (FTERM + DTERM)
  104. ENDIF
  105. C
  106. RETURN
  107. C------------- LAST LINE OF CHFIE FOLLOWS ------------------------------
  108. END