1 /* zzgfrelx.f -- translated by f2c (version 19980913).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include "f2c.h"
7
8 /* Table of constant values */
9
10 static integer c__5 = 5;
11 static integer c__7 = 7;
12 static integer c__0 = 0;
13 static logical c_true = TRUE_;
14
15 /* $Procedure ZZGFRELX ( Private --- GF, geometric relation finder ) */
zzgfrelx_(U_fp udstep,U_fp udrefn,U_fp udqdec,U_fp udcond,S_fp udfunc,char * relate,doublereal * refval,doublereal * tol,doublereal * adjust,doublereal * cnfine,integer * mw,integer * nw,doublereal * work,logical * rpt,S_fp udrepi,U_fp udrepu,S_fp udrepf,char * rptpre,char * rptsuf,logical * bail,L_fp udbail,doublereal * result,ftnlen relate_len,ftnlen rptpre_len,ftnlen rptsuf_len)16 /* Subroutine */ int zzgfrelx_(U_fp udstep, U_fp udrefn, U_fp udqdec, U_fp
17 udcond, S_fp udfunc, char *relate, doublereal *refval, doublereal *
18 tol, doublereal *adjust, doublereal *cnfine, integer *mw, integer *nw,
19 doublereal *work, logical *rpt, S_fp udrepi, U_fp udrepu, S_fp
20 udrepf, char *rptpre, char *rptsuf, logical *bail, L_fp udbail,
21 doublereal *result, ftnlen relate_len, ftnlen rptpre_len, ftnlen
22 rptsuf_len)
23 {
24 /* Initialized data */
25
26 static char cnames[80*7] = "< "
27 " " "= "
28 " " "> "
29 " "
30 " " "LOCMIN "
31 " " "ABSMIN "
32 " " "LOCMAX "
33 " "
34 " " "ABSMAX "
35 " ";
36 static logical cstep = FALSE_;
37
38 /* System generated locals */
39 integer work_dim1, work_dim2, work_offset, i__1, i__2, i__3;
40
41 /* Builtin functions */
42 integer s_cmp(char *, char *, ftnlen, ftnlen), s_rnge(char *, integer,
43 char *, integer);
44 /* Subroutine */ int s_copy(char *, char *, ftnlen, ftnlen);
45
46 /* Local variables */
47 doublereal addl, addr__;
48 integer case__;
49 logical need;
50 integer name__[2], pass, want;
51 doublereal step;
52 extern /* Subroutine */ int zzwninsd_(doublereal *, doublereal *, char *,
53 doublereal *, ftnlen);
54 integer i__;
55 extern integer cardd_(doublereal *);
56 extern /* Subroutine */ int zzgfwsts_(doublereal *, doublereal *, char *,
57 doublereal *, ftnlen), chkin_(char *, ftnlen), ucase_(char *,
58 char *, ftnlen, ftnlen), errch_(char *, char *, ftnlen, ftnlen);
59 integer minat;
60 doublereal endpt[2];
61 integer maxat;
62 doublereal value;
63 extern integer sized_(doublereal *);
64 extern /* Subroutine */ int errdp_(char *, doublereal *, ftnlen), copyd_(
65 doublereal *, doublereal *);
66 integer qcnum;
67 extern /* Subroutine */ int swapi_(integer *, integer *);
68 integer count;
69 doublereal start;
70 extern /* Subroutine */ int ljust_(char *, char *, ftnlen, ftnlen);
71 doublereal refer2;
72 extern logical failed_(void);
73 extern /* Subroutine */ int zzgfsolvx_(S_fp, U_fp, U_fp, U_fp, logical *,
74 L_fp, logical *, doublereal *, doublereal *, doublereal *,
75 doublereal *, logical *, U_fp, doublereal *), scardd_(integer *,
76 doublereal *);
77 extern integer isrchc_(char *, integer *, char *, ftnlen, ftnlen),
78 wncard_(doublereal *);
79 extern logical return_(void);
80 char contxt[500], locrel[80];
81 doublereal extrem, finish;
82 integer winsiz;
83 extern /* Subroutine */ int setmsg_(char *, ftnlen), errint_(char *,
84 integer *, ftnlen), sigerr_(char *, ftnlen), chkout_(char *,
85 ftnlen), ssized_(integer *, doublereal *), wnexpd_(doublereal *,
86 doublereal *, doublereal *), wnfetd_(doublereal *, integer *,
87 doublereal *, doublereal *), wnextd_(char *, doublereal *, ftnlen)
88 , wnintd_(doublereal *, doublereal *, doublereal *), wndifd_(
89 doublereal *, doublereal *, doublereal *), zzgfref_(doublereal *);
90
91 /* $ Abstract */
92
93 /* SPICE private routine intended solely for the support of SPICE */
94 /* routines. Users should not call this routine directly due to the */
95 /* volatile nature of this routine. */
96
97 /* This routine determines time intervals when the value of some */
98 /* geometric quantity related to one or more objects and an observer */
99 /* satisfies a user specified constraint within time intervals */
100 /* specified by the window CNFINE. */
101
102 /* Sister routine to ZZGFREL. Copy any edits to ZZGFREL or ZZGFRELX */
103 /* to the sister routine. */
104
105 /* $ Disclaimer */
106
107 /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
108 /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
109 /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
110 /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
111 /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
112 /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
113 /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
114 /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
115 /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
116 /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
117
118 /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
119 /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
120 /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
121 /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
122 /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
123 /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
124
125 /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
126 /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
127 /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
128 /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
129
130 /* $ Required_Reading */
131
132 /* SPK */
133 /* TIME */
134 /* NAIF_IDS */
135 /* FRAMES */
136
137 /* $ Keywords */
138
139 /* EPHEMERIS */
140 /* GEOMETRY */
141 /* SEARCH */
142
143 /* $ Declarations */
144 /* $ Abstract */
145
146 /* This file contains public, global parameter declarations */
147 /* for the SPICELIB Geometry Finder (GF) subsystem. */
148
149 /* $ Disclaimer */
150
151 /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
152 /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
153 /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
154 /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
155 /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
156 /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
157 /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
158 /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
159 /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
160 /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
161
162 /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
163 /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
164 /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
165 /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
166 /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
167 /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
168
169 /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
170 /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
171 /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
172 /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
173
174 /* $ Required_Reading */
175
176 /* GF */
177
178 /* $ Keywords */
179
180 /* GEOMETRY */
181 /* ROOT */
182
183 /* $ Restrictions */
184
185 /* None. */
186
187 /* $ Author_and_Institution */
188
189 /* N.J. Bachman (JPL) */
190 /* L.E. Elson (JPL) */
191 /* E.D. Wright (JPL) */
192
193 /* $ Literature_References */
194
195 /* None. */
196
197 /* $ Version */
198
199 /* - SPICELIB Version 2.0.0 29-NOV-2016 (NJB) */
200
201 /* Upgraded to support surfaces represented by DSKs. */
202
203 /* Bug fix: removed declaration of NVRMAX parameter. */
204
205 /* - SPICELIB Version 1.3.0, 01-OCT-2011 (NJB) */
206
207 /* Added NWILUM parameter. */
208
209 /* - SPICELIB Version 1.2.0, 14-SEP-2010 (EDW) */
210
211 /* Added NWPA parameter. */
212
213 /* - SPICELIB Version 1.1.0, 08-SEP-2009 (EDW) */
214
215 /* Added NWRR parameter. */
216 /* Added NWUDS parameter. */
217
218 /* - SPICELIB Version 1.0.0, 21-FEB-2009 (NJB) (LSE) (EDW) */
219
220 /* -& */
221
222 /* Root finding parameters: */
223
224 /* CNVTOL is the default convergence tolerance used by the */
225 /* high-level GF search API routines. This tolerance is */
226 /* used to terminate searches for binary state transitions: */
227 /* when the time at which a transition occurs is bracketed */
228 /* by two times that differ by no more than CNVTOL, the */
229 /* transition time is considered to have been found. */
230
231 /* Units are TDB seconds. */
232
233
234 /* NWMAX is the maximum number of windows allowed for user-defined */
235 /* workspace array. */
236
237 /* DOUBLE PRECISION WORK ( LBCELL : MW, NWMAX ) */
238
239 /* Currently no more than twelve windows are required; the three */
240 /* extra windows are spares. */
241
242 /* Callers of GFEVNT can include this file and use the parameter */
243 /* NWMAX to declare the second dimension of the workspace array */
244 /* if necessary. */
245
246
247 /* Callers of GFIDST should declare their workspace window */
248 /* count using NWDIST. */
249
250
251 /* Callers of GFSEP should declare their workspace window */
252 /* count using NWSEP. */
253
254
255 /* Callers of GFRR should declare their workspace window */
256 /* count using NWRR. */
257
258
259 /* Callers of GFUDS should declare their workspace window */
260 /* count using NWUDS. */
261
262
263 /* Callers of GFPA should declare their workspace window */
264 /* count using NWPA. */
265
266
267 /* Callers of GFILUM should declare their workspace window */
268 /* count using NWILUM. */
269
270
271 /* ADDWIN is a parameter used to expand each interval of the search */
272 /* (confinement) window by a small amount at both ends in order to */
273 /* accommodate searches using equality constraints. The loaded */
274 /* kernel files must accommodate these expanded time intervals. */
275
276
277 /* FRMNLN is a string length for frame names. */
278
279
280 /* FOVTLN -- maximum length for FOV string. */
281
282
283 /* Specify the character strings that are allowed in the */
284 /* specification of field of view shapes. */
285
286
287 /* Character strings that are allowed in the */
288 /* specification of occultation types: */
289
290
291 /* Occultation target shape specifications: */
292
293
294 /* Specify the number of supported occultation types and occultation */
295 /* type string length: */
296
297
298 /* Instrument field-of-view (FOV) parameters */
299
300 /* Maximum number of FOV boundary vectors: */
301
302
303 /* FOV shape parameters: */
304
305 /* circle */
306 /* ellipse */
307 /* polygon */
308 /* rectangle */
309
310
311 /* End of file gf.inc. */
312
313 /* $ Brief_I/O */
314
315 /* VARIABLE I/O DESCRIPTION */
316 /* -------- --- -------------------------------------------------- */
317 /* LBCELL P SPICELIB cell lower bound. */
318 /* NWREQ P Minimum number of workspace windows. */
319 /* UDSTEP I Name of the routine that computes and returns a */
320 /* time step. */
321 /* UDREFN I Name of the routine that computes a refined time. */
322 /* UDQDEC I Name of the routine that computes whether the */
323 /* scalar quantity is decreasing. */
324 /* UDCOND I Name of the routine that computes the scalar */
325 /* quantity condition with-respect-to the constraint. */
326 /* UDFUNC I The routine that computes the scalar quantity of */
327 /* interest. */
328 /* RELATE I Operator that either looks for an extreme value */
329 /* (max, min, local, absolute) or compares the */
330 /* scalar quantity value and a number. */
331 /* REFVAL I Value used as reference for scalar quantity */
332 /* condition. */
333 /* TOL I Convergence tolerance in seconds. */
334 /* ADJUST I Allowed variation for absolute extremal */
335 /* scalar conditions. */
336 /* CNFINE I Confinement window */
337 /* MW I Size of workspace windows. */
338 /* NW I Number of workspace windows. */
339 /* WORK I Array containing workspace windows */
340 /* RPT I Progress reporter on ( .TRUE.) or off ( .FALSE. ) */
341 /* UDREPI I Function that initializes progress reporting. */
342 /* UDREPU I Function that updates the progress report. */
343 /* UDREPF I Function that finalizes progress reporting. */
344 /* RPTPRE I Progress reporter beginning message. */
345 /* RPTSUF I Progress reporter ending message. */
346 /* BAIL I Logical indicating program interrupt monitoring. */
347 /* UDBAIL I Name of a routine that signals a program interrupt. */
348 /* RESULT I-O SPICE window containing results. */
349
350
351 /* $ Detailed_Input */
352
353 /* UDSTEP the routine that computes a time step in an attempt to */
354 /* find a transition of the scalar quantity. In the */
355 /* context of this routine's algorithm, a "transition" */
356 /* occurs where the scalar quantity value changes from */
357 /* "decreasing" to "not decreasing" or vice versa. */
358
359 /* This routine relies on UDSTEP returning step sizes */
360 /* small enough so that state transitions within the */
361 /* confinement window are not overlooked. There must */
362 /* never be two roots A and B separated by less than */
363 /* STEP, where STEP is the minimum step size returned by */
364 /* UDSTEP for any value of ET in the interval [A, B]. */
365
366 /* The calling sequence for UDSTEP is: */
367
368 /* CALL UDSTEP ( ET, STEP ) */
369
370 /* where: */
371
372 /* ET is the input start time from which the */
373 /* algorithm is to search forward for a state */
374 /* transition. ET is expressed as seconds past */
375 /* J2000 TDB. ET is a DOUBLE PRECISION number. */
376
377 /* STEP is the output step size. STEP indicates */
378 /* how far to advance ET so that ET and */
379 /* ET+STEP may bracket a state transition and */
380 /* definitely do not bracket more than one */
381 /* state transition. STEP is a DOUBLE */
382 /* PRECISION number. Units are TDB seconds. */
383
384 /* If a constant step size is desired, the routine */
385 /* GFSTEP may be used. This is the default option. */
386
387 /* UDREFN the routine that computes a refinement in the times */
388 /* that bracket a transition point. In other words, once */
389 /* a pair of times have been detected such that the system */
390 /* is in different states at each of the two times, UDREFN */
391 /* selects an intermediate time which should be closer to */
392 /* the transition state than one of the two known times. */
393 /* The calling sequence for UDREFN is: */
394
395 /* CALL UDREFN ( T1, T2, S1, S2, T ) */
396
397 /* where the inputs are: */
398
399 /* T1 a time when the system is in state S1. */
400
401 /* T2 a time when the system is in state S2. T2 */
402 /* is assumed to be larger than T1. */
403
404 /* S1 a logical indicating the state of the system */
405 /* at time T1. */
406
407 /* S2 a logical indicating the state of the system */
408 /* at time T2. */
409
410 /* UDREFN may use or ignore the S1 and S2 values. */
411
412 /* The output is: */
413
414 /* T a time to check for a state transition */
415 /* between T1 and T2. */
416
417 /* If a simple bisection method is desired, the routine */
418 /* GFREFN may be used. This is the default option. */
419
420 /* UDQDEC the routine that determines if the scalar quantity */
421 /* calculated by UDFUNC is decreasing. */
422
423 /* The calling sequence: */
424
425 /* CALL UDQDEC ( UDFUNC, ET, ISDECR ) */
426
427 /* where: */
428
429 /* ET a double precision value representing */
430 /* ephemeris time, expressed as seconds past */
431 /* J2000 TDB, at which to determine the time */
432 /* derivative of UDFUNC. */
433
434 /* ISDECR a logical return indicating whether */
435 /* or not the scalar value returned by UDFUNC */
436 /* is decreasing. ISDECR returns true if the */
437 /* time derivative of UDFUNC at ET is */
438 /* negative. */
439
440 /* UDCOND the routine that determines if UDFUNC satisfies */
441 /* some constraint condition at epoch ET. */
442
443 /* The calling sequence: */
444
445 /* CALL UDCOND ( UDFUNC, ET, IN_CON ) */
446
447 /* where: */
448
449 /* ET a double precision value representing */
450 /* ephemeris time, expressed as seconds past */
451 /* J2000 TDB, at which to evaluate UDFUNC. */
452
453 /* IN_CON a logical value indicating whether */
454 /* or not UDFUNC satisfies the constraint */
455 /* at ET (TRUE) or not (FALSE). */
456
457 /* UDFUNC the routine that returns the value of the scalar */
458 /* quantity of interest at time ET. The calling sequence */
459 /* for UDFUNC is: */
460
461 /* CALL UDFUNC ( ET, VALUE ) */
462
463 /* where: */
464
465 /* ET a double precision value representing */
466 /* ephemeris time, expressed as seconds past */
467 /* J2000 TDB, at which to determine the scalar */
468 /* value. */
469
470 /* VALUE the double precision value of the scalar */
471 /* quantity at ET. */
472
473 /* RELATE is a comparison operator, indicating the numeric */
474 /* constraint of interest. Values are: */
475
476 /* '>' value of scalar quantity greater than REFVAL. */
477
478 /* '=' value of scalar quantity equal to REFVAL. */
479
480 /* '<' value of scalar quantity less than REFVAL. */
481
482 /* ABSMAX-the scalar quantity is at an absolute */
483 /* maximum. */
484
485 /* ABSMIN-the scalar quantity is at an absolute */
486 /* minimum. */
487
488 /* LOCMAX-the scalar quantity is at an local maximum. */
489
490 /* LOCMIN-the scalar quantity is at an local minimum. */
491
492 /* REFVAL reference value for scalar quantity (in */
493 /* radians, radians/sec, km, or km/sec as appropriate). */
494
495 /* TOL is a tolerance value used to determine convergence of */
496 /* root-finding operations. TOL is measured in seconds */
497 /* and is greater than zero. */
498
499 /* ADJUST the amount by which the numerical quantity is */
500 /* allowed to vary from an absolute extremum. If ADJUST */
501 /* is non-zero, the resulting window contains */
502 /* intervals when the scalar quantity has */
503 /* values either between ABSMIN and ABSMIN + ADJUST */
504 /* or between ABSMAX and ABSMAX - ADJUST. ADJUST must */
505 /* not be negative. */
506
507 /* CNFINE is a SPICE window that confines the bounds of the */
508 /* search. Note that like all windows (see windows.req) */
509 /* CNFINE can contain multiple time intervals. See the */
510 /* Examples section for information on how to create this */
511 /* window. */
512
513 /* MW is the cell size of the windows in the workspace array */
514 /* WORK. */
515
516 /* NW is the number of windows in the workspace array WORK. */
517 /* NW must be at least as large as the parameter NWREQ. */
518
519 /* WORK is an array used to store workspace windows. This */
520 /* array has dimensions WORK (-5 : MW, NW). */
521
522 /* RPT is a logical variable which controls whether the */
523 /* progress reporter is on or off. The progress reporter */
524 /* writes to the user's terminal. */
525
526 /* UDREPI the routine that initializes a progress report. */
527 /* When progress reporting is enabled, UDREPI */
528 /* is called at the start of a search. The calling */
529 /* sequence of UDREPI is: */
530
531 /* UDREPI ( CNFINE, RPTPRE, RPTSUF ) */
532
533 /* DOUBLE PRECISION CNFINE ( LBCELL : * ) */
534 /* CHARACTER*(*) RPTPRE */
535 /* CHARACTER*(*) RPTSUF */
536
537 /* where */
538
539 /* CNFINE */
540
541 /* is the confinement window passed into ZZGFRELX, and */
542
543 /* RPTPRE */
544 /* RPTSUF */
545
546 /* are prefix and suffix strings used in the progress */
547 /* report: these strings are intended to bracket a */
548 /* representation of the fraction of work done. */
549
550 /* If the user has no progress reporting initialization */
551 /* routine, the SPICELIB routine GFREPI may be used. This */
552 /* is the default option. */
553
554 /* UDREPU the routine that updates the progress report for a */
555 /* search. The calling sequence of UDREPU is: */
556
557 /* UDREPU (IVBEG, IVEND, ET ) */
558
559 /* DOUBLE PRECISION ET */
560 /* DOUBLE PRECISION IVBEG */
561 /* DOUBLE PRECISION IVEND */
562
563 /* where ET is an epoch belonging to the confinement */
564 /* window, IVBEG and IVEND are the start and stop times, */
565 /* respectively of the current confinement window */
566 /* interval. The ratio of the measure of the portion */
567 /* of CNFINE that precedes ET to the measure of CNFINE */
568 /* would be a logical candidate for the search's */
569 /* completion percentage; however the method of */
570 /* measurement is up to the user. */
571
572 /* If the user has no progress reporting update routine, */
573 /* the SPICELIB routine GFREPU may be used. This is the */
574 /* default option. */
575
576 /* UDREPF the routine that finalizes a progress report. UDREPF */
577 /* has no arguments. */
578
579 /* If the user has no progress reporting finalizing */
580 /* routine, the SPICELIB routine GFREPF may be used. This */
581 /* is the default option. */
582
583 /* RPTPRE is an array of strings containing the prefixes of */
584 /* the output messages reported by the progress reporter. */
585 /* The Ith element of RPTPRE is the prefix for the */
586 /* message corresponding to the Ith traversal of the */
587 /* confinement window executed by this routine; such */
588 /* traversals are called "passes." The number of passes */
589 /* executed depends on the relational operator RELATE. */
590 /* Searches for local extrema and unadjusted absolute */
591 /* extrema require one pass; searches for adjusted */
592 /* absolute extrema, equalities, and inequalities require */
593 /* two passes. */
594
595 /* An example of the contents of RPTPRE for a distance */
596 /* equality search: */
597
598 /* RPTPRE(1) = 'Distance pass 1 of 2' */
599 /* RPTPRE(2) = 'Distance pass 2 of 2' */
600
601 /* RPTSUF is an array of strings containing the suffixes of */
602 /* the output messages reported by the progress reporter. */
603 /* The Ith element of RPTSUF is the suffix for the */
604 /* message corresponding to the Ith pass. */
605
606 /* An example of the contents of RPTSUF for a distance */
607 /* equality search: */
608
609 /* RPTSUF(1) = 'done.' */
610 /* RPTSUF(2) = 'done.' */
611
612 /* For this search, the complete progress report message */
613 /* for the Ith pass has the form */
614
615 /* 'Distance pass I of 2 xxx.xx% done.' */
616
617 /* BAIL is a logical indicating whether or not interrupt */
618 /* signaling is enabled. */
619
620 /* UDBAIL the routine that checks to see whether an interrupt */
621 /* signal has been issued from, e.g. the keyboard. If */
622 /* this capability is not to be used, a dummy function, */
623 /* GFBAIL must be supplied. */
624
625 /* RESULT is an initialized SPICE window. RESULT is large */
626 /* enough to hold all of the intervals, within the */
627 /* confinement window, on which the specified condition */
628 /* is met. */
629
630 /* $ Detailed_Output */
631
632 /* RESULT is a SPICE window containing the time intervals within */
633 /* the confinement window, over which the specified */
634 /* condition is met. */
635
636 /* RESULT is emptied before new values are assigned to */
637 /* it. */
638
639 /* $ Parameters */
640
641 /* LBCELL is the SPICELIB cell lower bound. */
642
643 /* NWREQ is the required number of workspace windows; the */
644 /* input argument NW must not be less than NWREQ. */
645
646 /* $ Exceptions */
647
648 /* 1) A negative value for ADJUST causes the routine to signal */
649 /* the error SPICE(VALUEOUTOFRANGE). A non-zero value for ADJUST */
650 /* when RELATE has any value other than "ABSMIN" or "ABSMAX", */
651 /* causes the routine to signal the error SPICE(INVALIDVALUE). */
652
653 /* 2) If an improper comparison operator is specified, the error */
654 /* SPICE(NOTRECOGNIZED) is signaled. */
655
656 /* 3) If TOL is not greater than zero, the error */
657 /* SPICE(VALUEOUTOFRANGE) will be signaled by routines called */
658 /* from this routine. */
659
660 /* 4) If the number of workspace windows is less than NWREQ, the */
661 /* error SPICE(TOOFEWWINDOWS) is signaled. */
662
663 /* 5) If the window size MW is less than 2, the error */
664 /* SPICE(INVALIDDIMENSION) will be signaled. */
665
666 /* 6) If the output SPICE window RESULT has insufficient capacity */
667 /* to contain the number of intervals on which the specified */
668 /* visibility condition is met, the error will be diagnosed */
669 /* by a routine in the call tree of this routine. If the result */
670 /* window has size less than 2, the error SPICE(WINDOWTOOSMALL) */
671 /* will be signaled by this routine. */
672
673 /* $ Files */
674
675 /* None. */
676
677 /* $ Particulars */
678
679 /* This routine determines time intervals when the value of some */
680 /* scalar quantity related to one or more objects and an observer */
681 /* satisfies a user specified constraint. It puts these times in a */
682 /* result window called RESULT. It does this by first finding */
683 /* windows when the quantity of interest is either */
684 /* monotonically increasing or decreasing. These windows are then */
685 /* manipulated to give the final result. Note that the determination */
686 /* of "=" involves finding intervals where the quantity is "less */
687 /* than" to a tolerance of TOL. This means that the end points of */
688 /* these intervals are within TOL of being equal to the value. */
689
690 /* $ Examples */
691
692 /* See GFEVNT. */
693
694 /* $ Restrictions */
695
696 /* The kernel files to be used by ZZGFRELX must be loaded (normally */
697 /* via the SPICELIB routine FURNSH) before ZZGFRELX is called. */
698
699 /* $ Literature_References */
700
701 /* None. */
702
703 /* $ Author_and_Institution */
704
705 /* N.J. Bachman (JPL) */
706 /* L.S. Elson (JPL) */
707 /* W.L. Taber (JPL) */
708 /* I.M. Underwood (JPL) */
709 /* E.D. Wright (JPL) */
710
711 /* $ Version */
712
713 /* - SPICELIB Version 1.1.1 08-DEC-2010 (EDW) */
714
715 /* Edit to replace term "schedule" with "window." Edit to */
716 /* Procedure line text. */
717
718 /* Argument list changed, removing the argument for */
719 /* the reference value update routine, UDQREF. Setting and */
720 /* updating the reference value now occurrs in this */
721 /* (ZZGFRELX) routine. */
722
723 /* Edit to Detailed I/O description for UDREPI, UDREPU, UDREPF, */
724 /* and UDBAIL. The descriptions stated the wrong name for the */
725 /* default GF functions corresponding to these arguments. */
726
727 /* - SPICELIB Version 1.1.0 16-FEB-2010 (EDW) */
728
729 /* Modified version of ZZGFREL. This version calls ZZGFSOLVX. */
730
731 /* - SPICELIB Version 1.0.0 21-FEB-2009 (NJB) (LSE) (WLT) (IMU) (EDW) */
732
733 /* -& */
734 /* $ Index_Entries */
735
736 /* numeric scalar quantity satisfies a condition */
737
738 /* -& */
739
740 /* SPICELIB functions */
741
742
743 /* Local parameters */
744
745
746 /* Workspace window indices: */
747
748
749 /* Number of supported comparison operators. */
750
751 /* One-letter alias for LBCELL to make references to the workspace */
752 /* array tolerable: */
753
754
755 /* Context string length: */
756
757
758 /* Local variables */
759
760
761 /* Saved variables */
762
763
764 /* Below we initialize the list of comparison operator names. */
765
766 /* Parameter adjustments */
767 work_dim1 = *mw + 6;
768 work_dim2 = *nw;
769 work_offset = work_dim1 - 5;
770
771 /* Function Body */
772
773 /* Set constant step parameter to .FALSE.. */
774
775
776 /* Standard SPICE error handling. */
777
778 if (return_()) {
779 return 0;
780 }
781 chkin_("ZZGFRELX", (ftnlen)8);
782
783 /* Make sure we have enough workspace windows. */
784
785 if (*nw < 5) {
786 setmsg_("The number of workspace windows (#) is less than the minimu"
787 "m #.", (ftnlen)63);
788 errint_("#", nw, (ftnlen)1);
789 errint_("#", &c__5, (ftnlen)1);
790 sigerr_("SPICE(TOOFEWWINDOWS)", (ftnlen)20);
791 chkout_("ZZGFRELX", (ftnlen)8);
792 return 0;
793 }
794
795 /* Make sure the workspace windows can contain at least one interval. */
796
797 if (*mw < 2) {
798 setmsg_("Workspace window size was #; size must be at least 2.", (
799 ftnlen)53);
800 errint_("#", mw, (ftnlen)1);
801 sigerr_("SPICE(INVALIDDIMENSION)", (ftnlen)23);
802 chkout_("ZZGFRELX", (ftnlen)8);
803 return 0;
804 }
805
806 /* Check the result window size. */
807
808 if (sized_(result) < 2) {
809 setmsg_("Result window size was #; size must be at least 2.", (ftnlen)
810 50);
811 i__1 = sized_(result);
812 errint_("#", &i__1, (ftnlen)1);
813 sigerr_("SPICE(INVALIDDIMENSION)", (ftnlen)23);
814 chkout_("ZZGFRELX", (ftnlen)8);
815 return 0;
816 }
817
818 /* Make sure the requested comparison is one we recognize. */
819
820 ljust_(relate, locrel, relate_len, (ftnlen)80);
821 ucase_(locrel, locrel, (ftnlen)80, (ftnlen)80);
822 qcnum = isrchc_(locrel, &c__7, cnames, (ftnlen)80, (ftnlen)80);
823 if (qcnum == 0) {
824 setmsg_("The comparison operator, # is not recognized. Supported qu"
825 "antities are: <, =, >, LOCMIN, ABSMIN, LOCMAX, ABSMAX.", (
826 ftnlen)113);
827 errch_("#", relate, (ftnlen)1, relate_len);
828 sigerr_("SPICE(NOTRECOGNIZED)", (ftnlen)20);
829 chkout_("ZZGFRELX", (ftnlen)8);
830 return 0;
831 }
832
833 /* Confirm ADJUST is non-negative. */
834
835 if (*adjust < 0.) {
836 setmsg_("ADJUST was #; must be non-negative.", (ftnlen)35);
837 errdp_("#", adjust, (ftnlen)1);
838 sigerr_("SPICE(VALUEOUTOFRANGE)", (ftnlen)22);
839 chkout_("ZZGFRELX", (ftnlen)8);
840 return 0;
841 }
842
843 /* Confirm ADJUST equals zero unless LOCREL (RELATE) has value */
844 /* "ABSMAX" or "ABSMIN." */
845
846 if (s_cmp(locrel, "ABSMIN", (ftnlen)80, (ftnlen)6) != 0 && s_cmp(locrel,
847 "ABSMAX", (ftnlen)80, (ftnlen)6) != 0) {
848 if (*adjust != 0.) {
849 setmsg_("ADJUST should have value zero for all comparison operat"
850 "ors except ABSMAX and ABSMIN", (ftnlen)83);
851 sigerr_("SPICE(INVALIDVALUE)", (ftnlen)19);
852 chkout_("ZZGFRELX", (ftnlen)8);
853 return 0;
854 }
855 }
856
857 /* If the confinement window is empty, the result window must */
858 /* be empty as well. In this case, there's not much to do. */
859
860 if (cardd_(cnfine) == 0) {
861 scardd_(&c__0, result);
862 chkout_("ZZGFRELX", (ftnlen)8);
863 return 0;
864 }
865
866 /* We need to set up several working windows, one each for */
867 /* increasing and decreasing windows, one for the confining */
868 /* window and one for copying. */
869
870 ssized_(mw, &work[(i__1 = (work_dim1 << 1) - 5 - work_offset) < work_dim1
871 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_"
872 , (ftnlen)768)]);
873 ssized_(mw, &work[(i__1 = work_dim1 - 5 - work_offset) < work_dim1 *
874 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_",
875 (ftnlen)769)]);
876 ssized_(mw, &work[(i__1 = work_dim1 * 3 - 5 - work_offset) < work_dim1 *
877 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_",
878 (ftnlen)770)]);
879 ssized_(mw, &work[(i__1 = (work_dim1 << 2) - 5 - work_offset) < work_dim1
880 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_"
881 , (ftnlen)771)]);
882 ssized_(mw, &work[(i__1 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
883 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_",
884 (ftnlen)772)]);
885 name__[0] = 2;
886 name__[1] = 1;
887 if (failed_()) {
888 chkout_("ZZGFRELX", (ftnlen)8);
889 return 0;
890 }
891
892 /* For equality constraints, we work with a somewhat expanded */
893 /* version of the confinement window so we can find equality */
894 /* solutions that lie on the boundary of the original confinement */
895 /* window. The expansion amount is ADDWIN. For other cases the */
896 /* expansion amount is set to zero. */
897
898 if (s_cmp(relate, "=", relate_len, (ftnlen)1) == 0) {
899 addl = .5;
900 addr__ = .5;
901 } else {
902 addl = 0.;
903 addr__ = 0.;
904 }
905 copyd_(cnfine, &work[(i__1 = work_dim1 * 3 - 5 - work_offset) < work_dim1
906 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_"
907 , (ftnlen)798)]);
908 wnexpd_(&addl, &addr__, &work[(i__1 = work_dim1 * 3 - 5 - work_offset) <
909 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1,
910 "zzgfrelx_", (ftnlen)799)]);
911 if (failed_()) {
912 chkout_("ZZGFRELX", (ftnlen)8);
913 return 0;
914 }
915
916 /* Set the reference value. */
917
918 zzgfref_(refval);
919
920 /* Make a local copy of the reference value. */
921
922 refer2 = *refval;
923
924 /* Set the pass number for progress reporting. */
925
926 pass = 1;
927
928 /* Initialize the work in progress reporter. */
929
930 if (*rpt) {
931 (*udrepi)(&work[(i__1 = work_dim1 * 3 - 5 - work_offset) < work_dim1 *
932 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfr"
933 "elx_", (ftnlen)826)], rptpre + (pass - 1) * rptpre_len,
934 rptsuf + (pass - 1) * rptsuf_len, rptpre_len, rptsuf_len);
935 }
936
937 /* Look up the size of the confinement window... */
938
939 count = wncard_(&work[(i__1 = work_dim1 * 3 - 5 - work_offset) <
940 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1,
941 "zzgfrelx_", (ftnlen)832)]);
942
943 /* Start the window that contains intervals when the quantity of */
944 /* interest is decreasing. The result will contain all intervals in */
945 /* (expanded) CNFINE when the selected scalar quantity function */
946 /* is decreasing, since this is how ZZGFSOLVX is configured. */
947
948 i__1 = count;
949 for (i__ = 1; i__ <= i__1; ++i__) {
950
951 /* Locate the bounds for the I'th interval of the confinement */
952 /* window. Results are accumulated in the WORK array. */
953
954 wnfetd_(&work[(i__2 = work_dim1 * 3 - 5 - work_offset) < work_dim1 *
955 work_dim2 && 0 <= i__2 ? i__2 : s_rnge("work", i__2, "zzgfre"
956 "lx_", (ftnlen)845)], &i__, &start, &finish);
957 zzgfsolvx_((S_fp)udfunc, (U_fp)udqdec, (U_fp)udstep, (U_fp)udrefn,
958 bail, (L_fp)udbail, &cstep, &step, &start, &finish, tol, rpt,
959 (U_fp)udrepu, &work[(i__2 = (work_dim1 << 1) - 5 -
960 work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 :
961 s_rnge("work", i__2, "zzgfrelx_", (ftnlen)847)]);
962 if (failed_()) {
963 chkout_("ZZGFRELX", (ftnlen)8);
964 return 0;
965 }
966 if (*bail) {
967 if ((*udbail)()) {
968 if (*rpt) {
969 (*udrepf)();
970 }
971 chkout_("ZZGFRELX", (ftnlen)8);
972 return 0;
973 }
974 }
975 }
976 if (*rpt) {
977 (*udrepf)();
978 }
979
980 /* Let's think about what we have now. We have the intervals in the */
981 /* confinement window when a value of some kind is decreasing. */
982
983 /* The left endpoints are points at which the quantity begins */
984 /* decreasing, thus they are times when the quantity is at a local */
985 /* maximum (at least in the interior of the confinement window). */
986
987 /* The right endpoints are where the quantity stops decreasing. Thus */
988 /* those endpoints in the interior of the confinement window are */
989 /* local minima of the quantity. */
990
991 /* The complement relative to the confinement window is the set of */
992 /* intervals within the confinement window for which the quantity is */
993 /* increasing. At the left endpoints of the complement the */
994 /* function is increasing. Thus the interior left endpoints are */
995 /* local minima within the confinement window. The interior right */
996 /* endpoints are local maxima within the confinement window. */
997
998 /* Moreover, to within our ability to detect local extrema, there */
999 /* are no local extrema within any of the intervals. Thus, the */
1000 /* function may be regarded as monotone within each of */
1001 /* the intervals of these windows. Thus for any desired value of the */
1002 /* quantity, there is at most one time within each of the intervals */
1003 /* that the desired value is achieved. */
1004
1005 if (s_cmp(locrel, "LOCMIN", (ftnlen)80, (ftnlen)6) == 0) {
1006
1007 /* We are interested in only interior minima of the quantity. */
1008 /* These occur at right endpoints of the intervals in TEMPW */
1009 /* that are interior points of CNFINE. First extract the right */
1010 /* endpoints. Then find those that are contained in the initial */
1011 /* confinement window, excluding endpoints. */
1012
1013 wnextd_("R", &work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1014 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1015 i__1, "zzgfrelx_", (ftnlen)913)], (ftnlen)1);
1016 zzgfwsts_(&work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1017 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1018 i__1, "zzgfrelx_", (ftnlen)915)], cnfine, "()", result, (
1019 ftnlen)2);
1020 chkout_("ZZGFRELX", (ftnlen)8);
1021 return 0;
1022 } else if (s_cmp(locrel, "LOCMAX", (ftnlen)80, (ftnlen)6) == 0) {
1023
1024 /* We are interested in only interior maxima of the quantity. */
1025 /* These occur at right endpoints of the intervals in TEMPW */
1026 /* that are interior points of CNFINE. */
1027
1028 wnextd_("L", &work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1029 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1030 i__1, "zzgfrelx_", (ftnlen)927)], (ftnlen)1);
1031 zzgfwsts_(&work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1032 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1033 i__1, "zzgfrelx_", (ftnlen)929)], cnfine, "()", result, (
1034 ftnlen)2);
1035 chkout_("ZZGFRELX", (ftnlen)8);
1036 return 0;
1037 }
1038
1039 /* We will need the intervals when the quantity of interest is */
1040 /* increasing in value. */
1041
1042 if (s_cmp(locrel, "ABSMIN", (ftnlen)80, (ftnlen)6) == 0 || s_cmp(locrel,
1043 "ABSMAX", (ftnlen)80, (ftnlen)6) == 0) {
1044
1045 /* We need an absolute max or min over the window CNFINE. */
1046 /* But we have decreasing values in WORK(B,DECRES). */
1047 /* Make a copy of WORK(B,DECRES) then compute the windows */
1048 /* of decreasing or increasing quantity over the window CNFINE. */
1049
1050 copyd_(&work[(i__1 = (work_dim1 << 1) - 5 - work_offset) < work_dim1 *
1051 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfr"
1052 "elx_", (ftnlen)947)], &work[(i__2 = (work_dim1 << 2) - 5 -
1053 work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 :
1054 s_rnge("work", i__2, "zzgfrelx_", (ftnlen)947)]);
1055 wnintd_(cnfine, &work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1056 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1057 i__1, "zzgfrelx_", (ftnlen)949)], &work[(i__2 = work_dim1 * 5
1058 - 5 - work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ?
1059 i__2 : s_rnge("work", i__2, "zzgfrelx_", (ftnlen)949)]);
1060 copyd_(&work[(i__1 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
1061 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfre"
1062 "lx_", (ftnlen)950)], &work[(i__2 = (work_dim1 << 1) - 5 -
1063 work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 :
1064 s_rnge("work", i__2, "zzgfrelx_", (ftnlen)950)]);
1065 wndifd_(cnfine, &work[(i__1 = (work_dim1 << 1) - 5 - work_offset) <
1066 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1067 i__1, "zzgfrelx_", (ftnlen)952)], &work[(i__2 = work_dim1 * 5
1068 - 5 - work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ?
1069 i__2 : s_rnge("work", i__2, "zzgfrelx_", (ftnlen)952)]);
1070 copyd_(&work[(i__1 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
1071 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfre"
1072 "lx_", (ftnlen)953)], &work[(i__2 = work_dim1 - 5 -
1073 work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 :
1074 s_rnge("work", i__2, "zzgfrelx_", (ftnlen)953)]);
1075
1076 /* Here's what we plan to do, we want to look over two windows */
1077 /* DECREASING and INCREASING to search for the absolute max or */
1078 /* min. We start with DECREASING. In this window the max is */
1079 /* always at the left endpoint, The min is at the right */
1080 /* endpoint. In the INCREASING window the min is at the LEFT */
1081 /* endpoint of an interval, the max is at the RIGHT endpoint of */
1082 /* an interval */
1083
1084 minat = 2;
1085 maxat = 1;
1086
1087 /* As yet we still need to compute our first extremum. */
1088
1089 need = TRUE_;
1090
1091 /* The extrema search is logically the same for both */
1092 /* maximum and minimum. We just need to keep track of */
1093 /* our extremum and when we find a more extreme value */
1094 /* replace it. DECREASING is first. */
1095
1096 for (case__ = 1; case__ <= 2; ++case__) {
1097 if (s_cmp(locrel, "ABSMIN", (ftnlen)80, (ftnlen)6) == 0) {
1098 want = minat;
1099 } else if (s_cmp(locrel, "ABSMAX", (ftnlen)80, (ftnlen)6) == 0) {
1100 want = maxat;
1101 }
1102 winsiz = wncard_(&work[(i__2 = name__[(i__1 = case__ - 1) < 2 &&
1103 0 <= i__1 ? i__1 : s_rnge("name", i__1, "zzgfrelx_", (
1104 ftnlen)991)] * work_dim1 - 5 - work_offset) < work_dim1 *
1105 work_dim2 && 0 <= i__2 ? i__2 : s_rnge("work", i__2,
1106 "zzgfrelx_", (ftnlen)991)]);
1107 i__1 = winsiz;
1108 for (i__ = 1; i__ <= i__1; ++i__) {
1109 wnfetd_(&work[(i__3 = name__[(i__2 = case__ - 1) < 2 && 0 <=
1110 i__2 ? i__2 : s_rnge("name", i__2, "zzgfrelx_", (
1111 ftnlen)995)] * work_dim1 - 5 - work_offset) <
1112 work_dim1 * work_dim2 && 0 <= i__3 ? i__3 : s_rnge(
1113 "work", i__3, "zzgfrelx_", (ftnlen)995)], &i__, endpt,
1114 &endpt[1]);
1115 (*udfunc)(&endpt[(i__2 = want - 1) < 2 && 0 <= i__2 ? i__2 :
1116 s_rnge("endpt", i__2, "zzgfrelx_", (ftnlen)998)], &
1117 value);
1118 if (failed_()) {
1119 chkout_("ZZGFRELX", (ftnlen)8);
1120 return 0;
1121 }
1122
1123 /* Initialize the extreme value. This step will */
1124 /* be executed on the first pass through the */
1125 /* DECREASING interval. */
1126
1127 if (need) {
1128 need = FALSE_;
1129 extrem = value;
1130 }
1131
1132 /* Check to see if current VALUE is more extreme than */
1133 /* EXTREM. */
1134
1135 if (s_cmp(locrel, "ABSMIN", (ftnlen)80, (ftnlen)6) == 0) {
1136 if (*adjust == 0. && value <= extrem) {
1137
1138 /* Let's save the epoch in case it's that of the */
1139 /* absolute min. Add this endpoint as a singleton */
1140 /* interval to the RESULT window. */
1141
1142 scardd_(&c__0, result);
1143 s_copy(contxt, "Saving current candidate epoch at wh"
1144 "ich an absolute minimum may occur.", (ftnlen)
1145 500, (ftnlen)70);
1146 zzwninsd_(&endpt[(i__2 = want - 1) < 2 && 0 <= i__2 ?
1147 i__2 : s_rnge("endpt", i__2, "zzgfrelx_", (
1148 ftnlen)1035)], &endpt[(i__3 = want - 1) < 2 &&
1149 0 <= i__3 ? i__3 : s_rnge("endpt", i__3,
1150 "zzgfrelx_", (ftnlen)1035)], contxt, result, (
1151 ftnlen)500);
1152 }
1153 extrem = min(extrem,value);
1154 } else {
1155 if (*adjust == 0. && value >= extrem) {
1156
1157 /* Let's save the epoch in case it's that of the */
1158 /* absolute max. Add this endpoint as a singleton */
1159 /* interval to the RESULT window. */
1160
1161 scardd_(&c__0, result);
1162 s_copy(contxt, "Saving current candidate epoch at wh"
1163 "ich an absolute maximum may occur.", (ftnlen)
1164 500, (ftnlen)70);
1165 zzwninsd_(&endpt[(i__2 = want - 1) < 2 && 0 <= i__2 ?
1166 i__2 : s_rnge("endpt", i__2, "zzgfrelx_", (
1167 ftnlen)1057)], &endpt[(i__3 = want - 1) < 2 &&
1168 0 <= i__3 ? i__3 : s_rnge("endpt", i__3,
1169 "zzgfrelx_", (ftnlen)1057)], contxt, result, (
1170 ftnlen)500);
1171 }
1172 extrem = max(extrem,value);
1173 }
1174 }
1175 if (failed_()) {
1176 chkout_("ZZGFRELX", (ftnlen)8);
1177 return 0;
1178 }
1179
1180 /* When we go to the next window, the min and max are at */
1181 /* opposite ends of the intervals. */
1182
1183 swapi_(&minat, &maxat);
1184 }
1185
1186 /* If the adjustment is zero, we're done. */
1187
1188 if (*adjust == 0.) {
1189 chkout_("ZZGFRELX", (ftnlen)8);
1190 return 0;
1191 }
1192
1193 /* We have a non-zero adjustment. we have the extreme value. Now */
1194 /* we need to find the epochs when the extreme value is achieved, */
1195 /* allowing for adjustment. */
1196
1197 if (s_cmp(locrel, "ABSMIN", (ftnlen)80, (ftnlen)6) == 0) {
1198 refer2 = extrem + *adjust;
1199 } else {
1200
1201 /* The only other possible value of LOCREL within this block */
1202 /* is 'ABSMAX'. */
1203
1204 refer2 = extrem - *adjust;
1205 }
1206
1207 /* If we reach this point, we need to re-establish the */
1208 /* original expanded coverage of 'DECREASING' and 'INCREASING'. */
1209
1210 copyd_(&work[(i__1 = (work_dim1 << 2) - 5 - work_offset) < work_dim1 *
1211 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfr"
1212 "elx_", (ftnlen)1114)], &work[(i__2 = (work_dim1 << 1) - 5 -
1213 work_offset) < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 :
1214 s_rnge("work", i__2, "zzgfrelx_", (ftnlen)1114)]);
1215 }
1216 wndifd_(&work[(i__1 = work_dim1 * 3 - 5 - work_offset) < work_dim1 *
1217 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfrelx_",
1218 (ftnlen)1118)], &work[(i__2 = (work_dim1 << 1) - 5 - work_offset)
1219 < work_dim1 * work_dim2 && 0 <= i__2 ? i__2 : s_rnge("work", i__2,
1220 "zzgfrelx_", (ftnlen)1118)], &work[(i__3 = work_dim1 - 5 -
1221 work_offset) < work_dim1 * work_dim2 && 0 <= i__3 ? i__3 : s_rnge(
1222 "work", i__3, "zzgfrelx_", (ftnlen)1118)]);
1223 if (failed_()) {
1224 chkout_("ZZGFRELX", (ftnlen)8);
1225 return 0;
1226 }
1227
1228 /* We have some kind of greater than, less than, or equal to */
1229 /* relation to solve for. Note that ABSMAX and ABSMIN are for case */
1230 /* where there is a non-zero adjustment. Reset the reference value, */
1231 /* which may have been changed in the ABSOLUTE MAX or MIN blocks */
1232 /* above. */
1233
1234 zzgfref_(&refer2);
1235
1236 /* If progress reporting is enabled, initialize the progress */
1237 /* reporter for a second pass over the confinement window. */
1238
1239 if (*rpt) {
1240
1241 /* Note that the window passed to UDREPI need not contain the */
1242 /* same intervals as those passed to UDREPU; the window passed to */
1243 /* UPREPI need only have the correct measure. From UDREPI's */
1244 /* perspective, the sole purpose of this window is to convey to */
1245 /* the progress reporting system the sum of the measures of the */
1246 /* increasing and decreasing windows. */
1247
1248 pass = 2;
1249 (*udrepi)(&work[(i__1 = work_dim1 * 3 - 5 - work_offset) < work_dim1 *
1250 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfr"
1251 "elx_", (ftnlen)1152)], rptpre + (pass - 1) * rptpre_len,
1252 rptsuf + (pass - 1) * rptsuf_len, rptpre_len, rptsuf_len);
1253 }
1254
1255 /* Find those intervals when the scalar quantity is less than */
1256 /* REFER2. */
1257
1258 scardd_(&c__0, result);
1259 for (case__ = 1; case__ <= 2; ++case__) {
1260 winsiz = wncard_(&work[(i__2 = name__[(i__1 = case__ - 1) < 2 && 0 <=
1261 i__1 ? i__1 : s_rnge("name", i__1, "zzgfrelx_", (ftnlen)1164)]
1262 * work_dim1 - 5 - work_offset) < work_dim1 * work_dim2 && 0
1263 <= i__2 ? i__2 : s_rnge("work", i__2, "zzgfrelx_", (ftnlen)
1264 1164)]);
1265
1266 /* Search each interval of the window identified by NAME(CASE) for */
1267 /* times when the quantity is less than the reference value. */
1268
1269 i__1 = winsiz;
1270 for (i__ = 1; i__ <= i__1; ++i__) {
1271 wnfetd_(&work[(i__3 = name__[(i__2 = case__ - 1) < 2 && 0 <= i__2
1272 ? i__2 : s_rnge("name", i__2, "zzgfrelx_", (ftnlen)1172)]
1273 * work_dim1 - 5 - work_offset) < work_dim1 * work_dim2 &&
1274 0 <= i__3 ? i__3 : s_rnge("work", i__3, "zzgfrelx_", (
1275 ftnlen)1172)], &i__, &start, &finish);
1276
1277 /* For each interval, accumulate the result in RESULT. */
1278
1279 /* Note we know that the behavior of the quantity is monotonic */
1280 /* within each window, so the step size can be large. In fact, */
1281 /* we use the interval length as the step size. */
1282
1283 step = finish - start;
1284 zzgfsolvx_((S_fp)udfunc, (U_fp)udcond, (U_fp)udstep, (U_fp)udrefn,
1285 bail, (L_fp)udbail, &c_true, &step, &start, &finish, tol,
1286 rpt, (U_fp)udrepu, result);
1287 if (failed_()) {
1288 chkout_("ZZGFRELX", (ftnlen)8);
1289 return 0;
1290 }
1291 if (*bail) {
1292 if ((*udbail)()) {
1293 chkout_("ZZGFRELX", (ftnlen)8);
1294 return 0;
1295 }
1296 }
1297 }
1298 }
1299 if (*rpt) {
1300
1301 /* Finish the progress report for the second pass. */
1302
1303 (*udrepf)();
1304 }
1305
1306 /* RESULT is the window, within the expanded confinement window, */
1307 /* over which the function of interest is less than the reference */
1308 /* value. We can use this window to get whatever was requested. */
1309
1310 if (s_cmp(locrel, "<", (ftnlen)80, (ftnlen)1) == 0 || s_cmp(locrel, "ABS"
1311 "MIN", (ftnlen)80, (ftnlen)6) == 0) {
1312
1313 /* We simply need to restrict our result to the original */
1314 /* confinement window. Note that the ABSMIN search with */
1315 /* non-zero adjustment is now a search for values less than the */
1316 /* adjusted absolute minimum. Same for ABSMAX below. */
1317
1318 wnintd_(cnfine, result, &work[(i__1 = work_dim1 * 5 - 5 - work_offset)
1319 < work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1320 i__1, "zzgfrelx_", (ftnlen)1224)]);
1321 copyd_(&work[(i__1 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
1322 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfre"
1323 "lx_", (ftnlen)1225)], result);
1324 } else if (s_cmp(locrel, ">", (ftnlen)80, (ftnlen)1) == 0 || s_cmp(locrel,
1325 "ABSMAX", (ftnlen)80, (ftnlen)6) == 0) {
1326
1327 /* Subtract from the confinement window the window where the */
1328 /* quantity is less than the reference value: the remainder is */
1329 /* the portion of the confinement window on which the quantity is */
1330 /* greater than or equal to the reference value. */
1331
1332 wndifd_(cnfine, result, &work[(i__1 = work_dim1 * 5 - 5 - work_offset)
1333 < work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1334 i__1, "zzgfrelx_", (ftnlen)1235)]);
1335 copyd_(&work[(i__1 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
1336 work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work", i__1, "zzgfre"
1337 "lx_", (ftnlen)1236)], result);
1338 } else {
1339
1340 /* This is the branch for the relational operator '='. */
1341
1342 /* Create a window of singleton intervals from the endpoints */
1343 /* of RESULT. */
1344
1345 scardd_(&c__0, &work[(i__1 = work_dim1 * 5 - 5 - work_offset) <
1346 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1347 i__1, "zzgfrelx_", (ftnlen)1245)]);
1348 i__1 = cardd_(result);
1349 for (i__ = 1; i__ <= i__1; ++i__) {
1350 s_copy(contxt, "Inserting endpoints of result window into worksp"
1351 "ace window WORK(B,TEMPW). These points are candidate epo"
1352 "chs that may satisfy an equality constraint.", (ftnlen)
1353 500, (ftnlen)148);
1354 zzwninsd_(&result[i__ + 5], &result[i__ + 5], contxt, &work[(i__2
1355 = work_dim1 * 5 - 5 - work_offset) < work_dim1 *
1356 work_dim2 && 0 <= i__2 ? i__2 : s_rnge("work", i__2,
1357 "zzgfrelx_", (ftnlen)1254)], (ftnlen)500);
1358 if (failed_()) {
1359 chkout_("ZZGFRELX", (ftnlen)8);
1360 return 0;
1361 }
1362 }
1363
1364 /* The window WORK(B,TEMPW) contains singleton intervals where */
1365 /* either the equality constraint is met, or where a boundary */
1366 /* point of the expanded confinement window is located. We're not */
1367 /* interested in the boundary points; these are likely not */
1368 /* solution points and in any case are outside the original */
1369 /* confinement window. */
1370
1371 /* Keep only the endpoints of RESULT that are contained in the */
1372 /* original confinement window CNFINE; these are by construction */
1373 /* interior points of the expanded confinement window. */
1374
1375 wnintd_(cnfine, &work[(i__1 = work_dim1 * 5 - 5 - work_offset) <
1376 work_dim1 * work_dim2 && 0 <= i__1 ? i__1 : s_rnge("work",
1377 i__1, "zzgfrelx_", (ftnlen)1276)], result);
1378 }
1379 chkout_("ZZGFRELX", (ftnlen)8);
1380 return 0;
1381 } /* zzgfrelx_ */
1382
1383