00001 C Work performed under the auspices of the U.S. Department of Energy 00002 C by Lawrence Livermore National Laboratory under contract number 00003 C W-7405-Eng-48. 00004 C 00005 SUBROUTINE DCNSTR (NEQ, Y, YNEW, ICNSTR, TAU, RLX, IRET, IVAR) 00006 C 00007 C***BEGIN PROLOGUE DCNSTR 00008 C***DATE WRITTEN 950808 (YYMMDD) 00009 C***REVISION DATE 950814 (YYMMDD) 00010 C 00011 C 00012 C----------------------------------------------------------------------- 00013 C***DESCRIPTION 00014 C 00015 C This subroutine checks for constraint violations in the proposed 00016 C new approximate solution YNEW. 00017 C If a constraint violation occurs, then a new step length, TAU, 00018 C is calculated, and this value is to be given to the linesearch routine 00019 C to calculate a new approximate solution YNEW. 00020 C 00021 C On entry: 00022 C 00023 C NEQ -- size of the nonlinear system, and the length of arrays 00024 C Y, YNEW and ICNSTR. 00025 C 00026 C Y -- real array containing the current approximate y. 00027 C 00028 C YNEW -- real array containing the new approximate y. 00029 C 00030 C ICNSTR -- INTEGER array of length NEQ containing flags indicating 00031 C which entries in YNEW are to be constrained. 00032 C if ICNSTR(I) = 2, then YNEW(I) must be .GT. 0, 00033 C if ICNSTR(I) = 1, then YNEW(I) must be .GE. 0, 00034 C if ICNSTR(I) = -1, then YNEW(I) must be .LE. 0, while 00035 C if ICNSTR(I) = -2, then YNEW(I) must be .LT. 0, while 00036 C if ICNSTR(I) = 0, then YNEW(I) is not constrained. 00037 C 00038 C RLX -- real scalar restricting update, if ICNSTR(I) = 2 or -2, 00039 C to ABS( (YNEW-Y)/Y ) < FAC2*RLX in component I. 00040 C 00041 C TAU -- the current size of the step length for the linesearch. 00042 C 00043 C On return 00044 C 00045 C TAU -- the adjusted size of the step length if a constraint 00046 C violation occurred (otherwise, it is unchanged). it is 00047 C the step length to give to the linesearch routine. 00048 C 00049 C IRET -- output flag. 00050 C IRET=0 means that YNEW satisfied all constraints. 00051 C IRET=1 means that YNEW failed to satisfy all the 00052 C constraints, and a new linesearch step 00053 C must be computed. 00054 C 00055 C IVAR -- index of variable causing constraint to be violated. 00056 C 00057 C----------------------------------------------------------------------- 00058 IMPLICIT DOUBLE PRECISION(A-H,O-Z) 00059 DIMENSION Y(NEQ), YNEW(NEQ), ICNSTR(NEQ) 00060 SAVE FAC, FAC2, ZERO 00061 DATA FAC /0.6D0/, FAC2 /0.9D0/, ZERO/0.0D0/ 00062 C----------------------------------------------------------------------- 00063 C Check constraints for proposed new step YNEW. If a constraint has 00064 C been violated, then calculate a new step length, TAU, to be 00065 C used in the linesearch routine. 00066 C----------------------------------------------------------------------- 00067 IRET = 0 00068 RDYMX = ZERO 00069 IVAR = 0 00070 DO 100 I = 1,NEQ 00071 C 00072 IF (ICNSTR(I) .EQ. 2) THEN 00073 RDY = ABS( (YNEW(I)-Y(I))/Y(I) ) 00074 IF (RDY .GT. RDYMX) THEN 00075 RDYMX = RDY 00076 IVAR = I 00077 ENDIF 00078 IF (YNEW(I) .LE. ZERO) THEN 00079 TAU = FAC*TAU 00080 IVAR = I 00081 IRET = 1 00082 RETURN 00083 ENDIF 00084 C 00085 ELSEIF (ICNSTR(I) .EQ. 1) THEN 00086 IF (YNEW(I) .LT. ZERO) THEN 00087 TAU = FAC*TAU 00088 IVAR = I 00089 IRET = 1 00090 RETURN 00091 ENDIF 00092 C 00093 ELSEIF (ICNSTR(I) .EQ. -1) THEN 00094 IF (YNEW(I) .GT. ZERO) THEN 00095 TAU = FAC*TAU 00096 IVAR = I 00097 IRET = 1 00098 RETURN 00099 ENDIF 00100 C 00101 ELSEIF (ICNSTR(I) .EQ. -2) THEN 00102 RDY = ABS( (YNEW(I)-Y(I))/Y(I) ) 00103 IF (RDY .GT. RDYMX) THEN 00104 RDYMX = RDY 00105 IVAR = I 00106 ENDIF 00107 IF (YNEW(I) .GE. ZERO) THEN 00108 TAU = FAC*TAU 00109 IVAR = I 00110 IRET = 1 00111 RETURN 00112 ENDIF 00113 C 00114 ENDIF 00115 100 CONTINUE 00116 00117 IF(RDYMX .GE. RLX) THEN 00118 TAU = FAC2*TAU*RLX/RDYMX 00119 IRET = 1 00120 ENDIF 00121 C 00122 RETURN 00123 C----------------------- END OF SUBROUTINE DCNSTR ---------------------- 00124 END
1.7.1
