text:
      SUBROUTINE LOOP_CT_CALLS_1(P,NHEL,H,IC)
C
C     Modules
C
      USE POLYNOMIAL_CONSTANTS
C
      IMPLICIT NONE
C
C     CONSTANTS
C
      INTEGER    NEXTERNAL
      PARAMETER (NEXTERNAL=4)
      INTEGER    NCOMB
      PARAMETER (NCOMB=16)
      INTEGER NBORNAMPS
      PARAMETER (NBORNAMPS=1)
      INTEGER    NLOOPS, NLOOPGROUPS, NCTAMPS
      PARAMETER (NLOOPS=35, NLOOPGROUPS=25, NCTAMPS=15)
      INTEGER    NLOOPAMPS
      PARAMETER (NLOOPAMPS=50)
      INTEGER    NWAVEFUNCS,NLOOPWAVEFUNCS
      PARAMETER (NWAVEFUNCS=6,NLOOPWAVEFUNCS=73)
      REAL*8     ZERO
      PARAMETER (ZERO=0D0)
      REAL*16     MP__ZERO
      PARAMETER (MP__ZERO=0.0E0_16)
C     These are constants related to the split orders
      INTEGER    NSO, NSQUAREDSO, NAMPSO
      PARAMETER (NSO=2, NSQUAREDSO=1, NAMPSO=2)
C
C     ARGUMENTS
C
      REAL*8 P(0:3,NEXTERNAL)
      INTEGER NHEL(NEXTERNAL), IC(NEXTERNAL)
      INTEGER H
C
C     LOCAL VARIABLES
C
      INTEGER I,J,K
      COMPLEX*16 COEFS(MAXLWFSIZE,0:VERTEXMAXCOEFS-1,MAXLWFSIZE)

      LOGICAL DUMMYFALSE
      DATA DUMMYFALSE/.FALSE./
C
C     GLOBAL VARIABLES
C
      INCLUDE 'coupl.inc'
      INCLUDE 'mp_coupl.inc'

      INTEGER HELOFFSET
      INTEGER GOODHEL(NCOMB)
      LOGICAL GOODAMP(NSQUAREDSO,NLOOPGROUPS)
      COMMON/FILTERS/GOODAMP,GOODHEL,HELOFFSET

      LOGICAL CHECKPHASE
      LOGICAL HELDOUBLECHECKED
      COMMON/INIT/CHECKPHASE, HELDOUBLECHECKED

      INTEGER SQSO_TARGET
      COMMON/SOCHOICE/SQSO_TARGET

      LOGICAL UVCT_REQ_SO_DONE,MP_UVCT_REQ_SO_DONE,CT_REQ_SO_DONE
     $ ,MP_CT_REQ_SO_DONE,LOOP_REQ_SO_DONE,MP_LOOP_REQ_SO_DONE
     $ ,CTCALL_REQ_SO_DONE,FILTER_SO
      COMMON/SO_REQS/UVCT_REQ_SO_DONE,MP_UVCT_REQ_SO_DONE
     $ ,CT_REQ_SO_DONE,MP_CT_REQ_SO_DONE,LOOP_REQ_SO_DONE
     $ ,MP_LOOP_REQ_SO_DONE,CTCALL_REQ_SO_DONE,FILTER_SO

      INTEGER I_SO
      COMMON/I_SO/I_SO
      INTEGER I_LIB
      COMMON/I_LIB/I_LIB

      COMPLEX*16 AMP(NBORNAMPS)
      COMMON/AMPS/AMP
      COMPLEX*16 W(20,NWAVEFUNCS)
      COMMON/W/W

      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
     $ -1:NLOOPWAVEFUNCS)
      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
      COMMON/WL/WL,PL

      COMPLEX*16 AMPL(3,NCTAMPS)
      COMMON/AMPL/AMPL

C
C     ----------
C     BEGIN CODE
C     ----------

C     The target squared split order contribution is already reached
C      if true.
      IF (FILTER_SO.AND.CTCALL_REQ_SO_DONE) THEN
        GOTO 1001
      ENDIF

C     CutTools call for loop numbers 1,29,30,2
      CALL LOOP_2(5,6,DCMPLX(CMASS_MDL_MW),DCMPLX(ZERO),2,I_SO,1)
C     CutTools call for loop numbers 3
      CALL LOOP_3(3,4,5,DCMPLX(ZERO),DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,2,I_SO,2)
C     CutTools call for loop numbers 4
      CALL LOOP_3(1,2,6,DCMPLX(ZERO),DCMPLX(ZERO),DCMPLX(ZERO),2,I_SO
     $ ,3)
C     CutTools call for loop numbers 5
      CALL LOOP_3(1,2,6,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW),DCMPLX(ZERO)
     $ ,2,I_SO,4)
C     CutTools call for loop numbers 6
      CALL LOOP_4(1,2,3,4,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,DCMPLX(ZERO),DCMPLX(ZERO),2,I_SO,5)
C     CutTools call for loop numbers 7
      CALL LOOP_3(1,2,6,DCMPLX(ZERO),DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,2,I_SO,6)
C     CutTools call for loop numbers 8
      CALL LOOP_4(1,2,4,3,DCMPLX(ZERO),DCMPLX(ZERO),DCMPLX(ZERO)
     $ ,DCMPLX(CMASS_MDL_MW),2,I_SO,7)
C     CutTools call for loop numbers 9,25
      CALL LOOP_1_2(2,5,6,DCMPLX(CMASS_MDL_MZ),0,I_SO,8)
C     CutTools call for loop numbers 10,26,31,32,11
      CALL LOOP_2(5,6,DCMPLX(CMASS_MDL_MW),DCMPLX(CMASS_MDL_MZ),2,I_SO
     $ ,9)
C     CutTools call for loop numbers 12
      CALL LOOP_3(3,4,5,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,DCMPLX(CMASS_MDL_MZ),2,I_SO,10)
C     CutTools call for loop numbers 13
      CALL LOOP_3(3,4,5,DCMPLX(CMASS_MDL_MZ),DCMPLX(ZERO),DCMPLX(ZERO)
     $ ,2,I_SO,11)
C     CutTools call for loop numbers 14
      CALL LOOP_3(3,4,5,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ)
     $ ,DCMPLX(CMASS_MDL_MW),2,I_SO,12)
C     CutTools call for loop numbers 15
      CALL LOOP_3(1,2,6,DCMPLX(CMASS_MDL_MZ),DCMPLX(ZERO),DCMPLX(ZERO)
     $ ,2,I_SO,13)
C     CutTools call for loop numbers 16
      CALL LOOP_3(1,2,6,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,DCMPLX(CMASS_MDL_MZ),2,I_SO,14)
C     CutTools call for loop numbers 17
      CALL LOOP_4(1,2,4,3,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ),2,I_SO,15)
C     CutTools call for loop numbers 18
      CALL LOOP_4(1,2,3,4,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW)
     $ ,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ),2,I_SO,16)
C     CutTools call for loop numbers 19
      CALL LOOP_3(1,2,6,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ)
     $ ,DCMPLX(CMASS_MDL_MW),2,I_SO,17)
C     CutTools call for loop numbers 20
      CALL LOOP_4(1,2,3,4,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ)
     $ ,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW),2,I_SO,18)
C     CutTools call for loop numbers 21
      CALL LOOP_4(1,2,4,3,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MZ)
     $ ,DCMPLX(ZERO),DCMPLX(CMASS_MDL_MW),2,I_SO,19)
C     CutTools call for loop numbers 22
      CALL LOOP_1_2(2,5,6,DCMPLX(CMASS_MDL_MH),0,I_SO,20)
C     CutTools call for loop numbers 23
      CALL LOOP_2(5,6,DCMPLX(CMASS_MDL_MW),DCMPLX(CMASS_MDL_MH),2,I_SO
     $ ,21)
C     CutTools call for loop numbers 24
      CALL LOOP_2(5,6,DCMPLX(CMASS_MDL_MH),DCMPLX(CMASS_MDL_MW),0,I_SO
     $ ,22)
C     CutTools call for loop numbers 27,28
      CALL LOOP_1_2(2,5,6,DCMPLX(CMASS_MDL_MW),0,I_SO,23)
C     CutTools call for loop numbers 33,35
      CALL LOOP_2(5,6,DCMPLX(ZERO),DCMPLX(ZERO),2,I_SO,24)
C     CutTools call for loop numbers 34
      CALL LOOP_2(5,6,DCMPLX(CMASS_MDL_MT),DCMPLX(ZERO),2,I_SO,25)
C     At this point, all reductions needed for (QCD=0 QED=6), i.e. of
C      split order ID=1, are computed.
      IF(FILTER_SO.AND.SQSO_TARGET.EQ.1) GOTO 5000

      GOTO 1001
 5000 CONTINUE
      CTCALL_REQ_SO_DONE=.TRUE.
 1001 CONTINUE
      END


text:
fl
\
fsclear
\
fswrite btstat.f

\ for the HC05
\ white wire
[ifndef hcProgPin
	25 wconstant hcProgPin
]
\ red wire
[ifndef hcRx
	26 wconstant hcRx
]
\  blue wire
[ifndef hcTx
	27 wconstant hcTx
]
\
[ifndef hcBaud
\ initital baud
\	d_9600 4/ wconstant hcBaud
\ operatingbaud
	d_230_400 4/ wconstant hcBaud
]
[ifndef hcSerialCog
	0 wconstant hcSerialCog
]
\
[ifndef hcSend
: hcSend
	cogid hcSerialCog iolink .cstr 1000 delms cogid iounlink cr
;
]
\
hcSerialCog cogreset 100 delms
c" hcRx hcTx hcBaud serial" hcSerialCog cogx 100 delms
1 hcSerialCog sersetflags
\
\
: hcSend
	cogid hcSerialCog iolink .cstr 1000 delms cogid iounlink cr
;
\
hcProgPin pinhi hcProgPin pinout
\
c" AT+VERSION~h0D~h0A" hcSend
c" AT+ROLE?~h0D~h0A" hcSend
c" AT+UART?~h0D~h0A" hcSend
c" AT+NAME?~h0D~h0A" hcSend
c" AT+PSWD?~h0D~h0A" hcSend
\
hcProgPin pinlo
\

...

fswrite btinit.f

\ for the HC05
\ white wire
[ifndef hcProgPin
	25 wconstant hcProgPin
]
\ red wire
[ifndef hcRx
	26 wconstant hcRx
]
\  blue wire
[ifndef hcTx
	27 wconstant hcTx
]

[ifndef hcBaud
\ initital baud
\	d_9600 4/ wconstant hcBaud
\ operatingbaud
	d_230_400 4/ wconstant hcBaud
]
[ifndef hcSerialCog
	0 wconstant hcSerialCog
]
\
[ifndef hcSend
: hcSend
	cogid hcSerialCog iolink .cstr 1000 delms cogid iounlink cr
;
]
\
hcSerialCog cogreset 100 delms
c" hcRx hcTx hcBaud serial" hcSerialCog cogx 100 delms
1 hcSerialCog sersetflags
\
: hcSend
	cogid hcSerialCog iolink .cstr 1000 delms cogid iounlink cr
;
\
hcProgPin pinhi hcProgPin pinout
\
c" AT+VERSION~h0D~h0A" hcSend
c" AT+ROLE?~h0D~h0A" hcSend
c" AT+UART?~h0D~h0A" hcSend
c" AT+NAME?~h0D~h0A" hcSend
c" AT+PSWD?~h0D~h0A" hcSend
c" AT+NAME=SIGHT_848~h0D~h0A" hcSend
c" AT+PSWD=8488~h0D~h0A" hcSend
c" AT+UART=230400,0,0~h0D~h0A" hcSend
c" AT+VERSION~h0D~h0A" hcSend
c" AT+ROLE?~h0D~h0A" hcSend
c" AT+UART?~h0D~h0A" hcSend
c" AT+NAME?~h0D~h0A" hcSend
c" AT+PSWD?~h0D~h0A" hcSend
\
hcProgPin pinlo
\

...

\
fswrite boot.f
hA state orC! version W@ .cstr cr cr cr
\ findEETOP ( -- n1 ) the top of the eeprom + 1
: findEETOP
\
\ search 32k block increments until we get a fail
\
	0
	h100000 h8000
	do
		i t0 2 eereadpage
		if
			leave
		else
			i h7FFE + t0 3 eereadpage
			if
				leave
			else
				drop i h8000 +
			then
		then
	h8000 +loop
;
\
c" boot.f - Finding top of eeprom, " .cstr findEETOP ' fstop 2+ alignl L! forget findEETOP c" Top of eeprom at: " .cstr fstop . cr
\
\
c" boot.f - DONE PropForth Loaded~h0D~h0D" .cstr hA state andnC!
\
\
1 wconstant build_servo
\
\
\ waitcnt ( n1 n2 -- n1 ) \ wait until n1, add n2 to n1
[ifndef waitcnt
: waitcnt
	_xasm2>1 h1F1 _cnip
;
]
\
\ a cog special register
[ifndef ctra
h1F8	wconstant ctra
]
\
\ a cog special register
[ifndef ctrb
h1F9	wconstant ctrb
]
\
\ a cog special register
[ifndef frqa
h1FA	wconstant frqa
]
\
\ a cog special register
[ifndef frqb
h1FB	wconstant frqb
]
\
\ a cog special register
[ifndef phsa
h1FC	wconstant phsa
]
\
\ a cog special register
[ifndef phsb
h1FD	wconstant phsb
]
\
\
\ Servo motors are controlled by a pulse width modulation (pwm).
\ Every 20 milliseconds a pulse of length 0.75 milliseconds
\ to 2.25 millseconds determines the position of the servo motor.
\
\ Example1: 0.75ms
\ ---_____________________________________________________________________________---________
\
\ Example2: 1.5ms
\ ------__________________________________________________________________________------_____
\
\ Example3: 2.25ms
\ ---------_______________________________________________________________________---------__
\
\
\ Each servo motor requires calibration, so we will set the absolute minimum to 0.5 ms and the
\ absolute maximum to 2.5 ms. Depending on the servos you are using, you may want to adjust these
\ values.
\

\
\ This driver will drive 16 servos per cog, it uses the 2 counters in a time domain multiplexed
\ mode to generate the pulses.
\
\ Each counter is set to single ended PWM/NCO mode. In this mode the counter drives an io pin.
\ The value of the pin is controlled by bit 31 of the phsa/phsb registers. These registers and
\ incremented every cycle of the system clock. For an 80 Mhz system clock this provides for
\ 12.5 ns resolution.
\
\ Without intervention the phsa/phsb registers would generate a square wave with the cycle time
\ being more than 50 seconds
\ ----------------------____________________----------------------____________________
\
\ By starting the phsa/phsb register at the appropriate place in the cycle, and by resetting
\ it to the same position every 20 ms, we generate a pwm signal on the io pin.
\
\ Since there are 2 counters we can do this with 2 pins at a time.
\
\ And since we only need intervention during the hi time of the signal, we can interleave the
\ processing so that each clock will address 8 io pins.
\
\ This is done by taking the 20ms cycle, and splitting it into 8 2.5ms sections, and driving the
\ corresponding pin hi in that portion of the cycle.
\
\ This means the absolute maximum we can set the pulse width will be 2.5ms
\
\ 0ms                 5ms                 10ms                15ms                20ms
\ ------__________________________________________________________________________------_____ ctra pin n1
\ __________------___________________________________________________________________________ ctra pin n1 + 2
\ ____________________------_________________________________________________________________ ctra pin n1 + 4
\ ______________________________------_______________________________________________________ ctra pin n1 + 6
\ ________________________________________------_____________________________________________ ctra pin n1 + 8
\ __________________________________________________------___________________________________ ctra pin n1 + 10
\ ____________________________________________________________------_________________________ ctra pin n1 + 12
\ ______________________________________________________________________------_______________ ctra pin n1 + 14
\
\ ------__________________________________________________________________________------_____ ctrb pin n1 + 1
\ __________------___________________________________________________________________________ ctrb pin n1 + 3
\ ____________________------_________________________________________________________________ ctrb pin n1 + 5
\ ______________________________------_______________________________________________________ ctrb pin n1 + 7
\ ________________________________________------_____________________________________________ ctrb pin n1 + 9
\ __________________________________________________------___________________________________ ctrb pin n1 + 11
\ ____________________________________________________________------_________________________ ctrb pin n1 + 13
\ ______________________________________________________________________------_______________ ctrb pin n1 + 15
\
\
\
\
\ This constant defines the number of times per second a pulse is generated
\ for the servo motor.
\
d_50 wconstant sm_cyclefreq
\
\
\ The number of clock cycles in one pulse sequence of 20ms
\
clkfreq sm_cyclefreq u/ constant sm_cyclecnt
\
\
\ The number of cycles in 2.5ms, 1/8 of the 20 ms cycle. However an effect
\ of this is that the maximum pulse time is 2.5ms less the time it take to
\ process a cycle. Thus the absolute maximum time it for a pulse width is
\ approximately 2.39ms to 2.4ms on an 80mHz system.
\
sm_cyclecnt 8 u/ constant sm_cyclecnt/8
\
\
\ These constants define the absolute minium and maximum lengths of time
\ for the pulse. The maximum can never be more than the time slot.
\
\ The minimum was set to 500usec and and maximum to 2500usec
\ but someone managed to turn their servo into a countinuos rotation servo
\ without intending to do so.
\ So the defaults are set to 750usec and 2250usec.
\
clkfreq d_750   d_1_000_000 u*/ constant sm_minhi
clkfreq d_2_250 d_1_000_000 u*/ sm_cyclecnt/8 min constant sm_maxhi
\
\
\ This variable defines which IO pins are enabled as outputs
\ by the servo driver, the default is 0 - 27 since 30 & 31 are used
\ for the serial port, and 28 & 29 for the eeprom
\
\ variable sm_enable h_0FFF_FFFF sm_enable L!
variable sm_enable h_0000_0060 sm_enable L!
\
\
\ this array defines the pulse hi times, it is initialized 1/2 way between
\ the min and max
\
lockdict variable sm_hitime d_124 allot freedict		\ 32 longs
\
\
\ this array defines the minimum pulse hi time, it is initialized to the minimum
\
lockdict variable sm_minhitime d_124 allot freedict	\ 32 longs
\
\
\ this array defines the maximum pulse hi time, it is initialized to the maximum
\
lockdict variable sm_maxhitime d_124 allot freedict	\ 32 longs
\
\
\ _sm_idx ( n1 addr -- addr ) calculates the array offset for the particular servo
: _sm_idx
	swap 0 max d_31 min 4* +
;
\
\
\ sm_setpos ( n1 u -- ) n1 is an integer between 0 and 31, u is an unsigned integer between 0 and 10,000
: sm_setpos
	0 max d_10_000 min
\						\ ( n1 u -- )
	over dup sm_minhitime _sm_idx L@
\						\ ( n1 u n1 minhitime -- )
	swap sm_maxhitime _sm_idx L@
\						\ ( n1 u minhitime maxhitime -- )
	over -
\						\ ( n1 u minhitime range -- )
	rot d_10_000 u*/ +
\						\ ( n1 pos -- )
	swap sm_hitime _sm_idx L!
;
\
\
\
\
\ sm_servo ( n1 -- ) this cog will drive servos n1 to n1 + 15
: sm_servo
	4 state andnC!
	c" SERVO" cds W!
\				\ ( n1 -- )
	0 max h18 min
	dup h_10 bounds
	do
		sm_enable L@ 1 i lshift and
		if
			i dup pinlo pinout
		then
	loop
\				\ ( n1 -- )
	1 frqa COG!
	1 frqb COG!
\				\ counter set to drive pin n1, single ended nco/pwm mode
	dup h_1000_0000 +
\				\ the offset into the array defining the pulse width
	swap 4* sm_hitime +
\				\ ( ctrn1 hitimeoffset -- )
	cnt COG@ sm_cyclecnt +
\				\ ( ctrn1 hitimeoffset nextcycletime -- )
	begin
		h10 0
		do
			0 phsa COG! 0 phsb COG!
\				\ counter a set to drive pin n1,   single ended nco/pwm mode
\				\ counter b set to drive pin n1+1, single ended nco/pwm mode
\				\ ( ctrn1 hitimeoffset nextcycletime -- )
			rot dup i + dup ctra COG!
			1 + ctrb COG!
\				\ ( hitimeoffset nextcycletime ctrn1 -- )
			rot dup i 2 lshift + dup L@ negate phsa COG!
			4 + L@ negate phsb COG!
\				\ ( nextcycletime ctrn1 hitimeoffset -- )
			rot
\				\ ( hitimeoffset ctrn1 nextcycletime -- )
\				\ wait for the next 2.5ms time slot
			sm_cyclecnt/8 waitcnt
		2 +loop
	0 until
;
\
\ _pos? ( n1 n2 -- n3 n4 ) n3 usec per cycle for pin n1, n4 usec per cycle for pin n2
: _pos?
\					\ increment counter when pulse is hi for pin n1
	h_6800_0000 + ctra COG! 1 frqa COG!
\					\ increment counter when pulse is hi for pin n2
	h_6800_0000 + ctrb COG! 1 frqb COG!
\					\ zero counts, wait one cycle and get counts
	cnt COG@ 0 phsb COG! 0 phsa COG!
	sm_cyclecnt + 0 waitcnt
	phsb COG@ phsa COG@
\					\ display number of usec pin was hi in one cycle count
	rot drop
\					( n1 n2 timen1 timen2 -- )
;

\ pos? ( n1 n2 -- ) display for pin n1 and pin n2 the number of usec hi per cycle
: pos?
	2dup
	_pos?
\
	>r
\					( n1 n2 timen1 -- )
 	rot  <# # # #> .cstr h2D emit hF4240 clkfreq u*/ <# # # # # #> .cstr space
	r>
 	swap <# # # #> .cstr h2D emit hF4240 clkfreq u*/ <# # # # # #> .cstr space
;
\
\
\ sm_setminmax ( min max n1 -- ) sets calibration parameters for each servo
: sm_setminmax
	swap sm_maxhi min over sm_maxhitime _sm_idx L!
	swap sm_minhi max swap sm_minhitime _sm_idx L!
;
\
\ _sm_cal2 ( n1 -- )
\ _sm_cal2
: _sm_cal2
	sm_enable L@ over 3 swap lshift and
	if
		dup 1+ pos?
		hD emit
	else
		drop
	then
;
\
\ _sm_cal1 ( servo key delta -- servo key )
: _sm_cal1
	rot tuck
\					\ ( key servo delta servo -- )
	sm_hitime _sm_idx dup
\					\ ( key servo delta addr addr -- )
	L@ rot + sm_minhi max sm_maxhi min swap L!
\					\ ( key servo -- )
	swap
;
\
\ _sm_cal ( n1 -- n1 ) n1 - the serv0 0 - 9
: _sm_cal
	begin
		key
		dup h61 = if h-2710 _sm_cal1 else
		dup h73 = if h-3E8 _sm_cal1 else
		dup h64 = if h-64 _sm_cal1 else
		dup h66 = if h-10 _sm_cal1 else
		dup h68 = if h10 _sm_cal1 else
		dup h6A = if h64 _sm_cal1 else
		dup h6B = if h3E8 _sm_cal1 else
		dup h6C = if h2710 _sm_cal1
		thens
		over _sm_cal2
		hD =
	until
;
\
\
\ sm_calibrate ( n1 -- ) n1 - servo
: sm_calibrate
	begin
		." a - <<<<    s - <<<    d - <<    f - <    h - >    j - >>    k - >>>    l - >>>>" cr

		." Move servo to leftmost position then hit enter" cr
		_sm_cal cr
		dup sm_hitime _sm_idx L@ swap
\										\ ( left n1 -- )
\
		." Move servo to rightmost position then hit enter" cr
		_sm_cal cr
		dup sm_hitime _sm_idx L@ swap
\										\ ( left right n1 -- )
		." Calibration parameters: " rot dup . rot dup . rot dup . cr
		dup >r sm_setminmax
		r>
\										\ ( center left-center right n1 -- )
\
		." 0 - Leftmost   1 - Center    2 - Rightmost     r - Recalibrate    ESC - Done" cr
		0 begin
			drop
			key
			dup h30 = if over 0 sm_setpos else
			dup h31 = if over h1388 sm_setpos else
			dup h32 = if over h2710 sm_setpos
			thens
			over _sm_cal2
\
			dup h1B = over 72 = or
		until
		h1B =
	until
	cr dup sm_minhitime _sm_idx L@ . dup sm_maxhitime _sm_idx L@ . . cr
;
\
\
\
\
\ sm_start_servos ( -- ) runs 32 servo drivers on cogs 0 - 1, io pins 0 - 31
\ 			modify after calibrating servos
: sm_start_servos
\ initialize arrays, default calibration
\
	32 0
	do
		sm_minhi i sm_minhitime _sm_idx L!
		sm_maxhi i sm_maxhitime _sm_idx L!
	loop
\
\ calibration for servos
\	105000 200000 0 sm_setminmax

	32 0
	do
		i sm_minhitime _sm_idx L@
		i sm_maxhitime _sm_idx L@
		over - 2/ +
\
		i sm_hitime _sm_idx L!
	loop
\
	1 cogreset h_10 delms c" 0 sm_servo" 1 cogx
\	2 cogreset h_10 delms c" d_16  sm_servo" 2 cogx
	10 delms
;
\
\
\
1 wconstant build_sightbox
: sbv c" Sight Box Software Version 2013-JAN-15 11:58:00~h0D~h0D" ;
\
\
\
d_40 wconstant cal_time
\
fsload cal.f
\
\
\ Set the console bluetooth serial port
\
\
\
\ running pin definitions
\
\
\ bluetooth
\
\ white wire
25 wconstant hcProgPin
\
\ red wire
26 wconstant hcRx
\
\  blue wire
27 wconstant hcTx
\
\ led on board
0 wconstant _red_led
1 wconstant _green_led
2 wconstant _yellow_led
\
4 wconstant _laser
5 wconstant _vertical_servo
6 wconstant _horizontal_servo
\
\ 8 pin connectoreon right - from rear forward
\
\ 5v   - red
\ 3.3v - nc
\ GND  - brown
\ IO-7 - nc
\ IO-6 - yellow
\ IO-5 - green
\ IO-4 - white
\ IO-3 - nc
\
\
\
\
\ RJ45 connector
\
\
\ Orange/White - GND
\ Orange       - IO-09 -- Pullup / remote red led
\ Green/White  - IO-10 -- Top switch PB
\ Green        - IO-11 -- Top switch rotate 1
\ Brown/White  - IO-12 -- Top switch rotate 2
\ Brown        - IO-13 -- Side switch PB
\ Blue/White   - IO-14 -- Side switch rotate 1
\ Blue         - IO-15 -- Side switch rotate 2
\
\
\ Bluetooth
\
\ ID: SIGHT_848
\ Pairing Code: 8488
\
\ Debug bits
\
\ bit 0 - Echo rotate and up down commands
\ bit 1 - Echo tick command
\ bit 2 - Display in port bit changes
\ bit 3 - Display XY commands
\ bit 4 - Debug information - display every 8 seconds
\
\
wvariable in_debug 0 in_debug W!
variable bootTicks
variable _in_nextTickCount
\
\ _in_monitor1 ( 0 old new diff -- cmd old new diff -- )
: _in_monitor1
	dup h2000 and if
		over h2000 and
		if
			c" sideUp" 3 ST!
		else
			c" sideDown" 3 ST!
		then
	else dup h8000 and if
		over h8000 and 0=
		if
			over h4000 and
			if
				c" sideCW" 3 ST!
			else
				c" sideCCW" 3 ST!
			then
		then
	else dup h0400 and if
		over h0400 and
		if
			c" topUp" 3 ST!
		else
			c" topDown" 3 ST!
		then
	else dup h1000 and if
		over h1000 and 0=
		if
			over h0800 and
			if
				c" topCW"  3 ST!
			else
				c" topCCW"  3 ST!
			then
		then
	thens
;
\
\ _in_monitor2 ( cmd -- )
: _in_monitor2
	dup if
		in_debug W@ 1 and
		if
			c" DEBUG: " .concstr dup .concstr .concr
		then
		.cstr cr
	else
		drop
	then
;
\
: in_monitor
	4 state andnC!
	c" IN_MONITOR" cds W!
	9 pinout 9 pinhi ina COG@ hFC00 and
	cnt COG@ _in_nextTickCount L!
	0 bootTicks L!
	begin

		bootTicks L@ cal_time > cnt COG@ or h_400_0000 and 9 px
\ (old -- )
		9 px?
		if
			ina COG@ hFC00 and
			2dup xor 0<>
			if
\ ( old new -- )

				0 rot2 2dup xor
\ ( cmd old new diff -- )
				in_debug W@ 4 and
				if
					base W@ 2 base W!
					3 ST@ .conword  bl .conemit
					2 ST@ .conword  bl .conemit
					over .conword .concr
					base W!
				then

				_in_monitor1
\ ( cmd old new diff -- )
				drop rot
				_in_monitor2
			then
\ ( old new -- )
			nip
		then
\ ( new -- )
		_in_nextTickCount L@ cnt COG@ - 0 <=
		if
			c" tick"
			in_debug W@ 2 and
			if
				c" DEBUG: " .concstr dup .concstr .concr
			then
			.cstr cr

			_in_nextTickCount L@ clkfreq 4/ + _in_nextTickCount L!
		then
		0
	until
	drop
;
\
wvariable in_topDown
wvariable in_sideDown
wvariable in_laserCount
wvariable in_scan
wvariable in_calibrate
wvariable in_laserFlash
\
wvariable x
wvariable y
wvariable minX
wvariable minY
wvariable maxX
wvariable maxY
\
wvariable _clickIncrement
\
wvariable _scaledxy
\
\
d_333 wconstant _defClickIncrement
\
\
: _dbg
	c" ~h0D~h0Din_topDown: " .concstr in_topDown W@ .con
	c" ~h0Din_sideDown: " .concstr in_sideDown W@ .con
	c" ~h0Din_laserCount: " .concstr in_laserCount W@ .con
	c" ~h0Din_scan: " .concstr in_scan W@ .con
	c" ~h0Din_calibrate: " .concstr in_calibrate W@ .con
	c" ~h0Din_laserFlash: " .concstr in_laserFlash W@ .con
	c" ~h0Dx: " .concstr x W@ .con
	c" ~h0Dy: " .concstr y W@ .con
	c" ~h0DminX: " .concstr minX W@ .con
	c" ~h0DminY: " .concstr minY W@ .con
	c" ~h0DmaxX: " .concstr maxX W@ .con
	c" ~h0DmaxY: " .concstr maxY W@ .con
	c" ~h0D_clickIncrement: " .concstr _clickIncrement W@ .con
	c" ~h0D_scaledxy: " .concstr _scaledxy W@ .con
	c" ~h0DbootTicks: " .concstr bootTicks L@ .con
	c" ~h0D_in_nextTickCount: " .concstr _in_nextTickCount L@ .con
	c" ~h0D_in_debug: " .concstr in_debug L@ .con
	c" ~h0Dcal_time: " .concstr  cal_time .con .concr .concr
;
\ limit ( x -- X )
: limit
	d_9999 min 0 max
;
\
\ _scalexy( x y -- scaleX scaleY)
: _scalexy
	limit swap limit swap
	d_2000 /mod
\ ( x ymod ydiv -- )
	rot d_2000 /mod
\ ( ymod ydiv xmod xdiv -- )
	4* cal_table +
\ ( ymod ydiv xmod caloff -- )
	rot d_24 u* +
\ ( ymod xmod caloff -- )
\	over . dup W@ . dup 4+ W@ . cr
\	2 ST@ . dup 2+ W@ . dup d_26 + W@ . cr

	dup 4+ W@ over W@ - rot d_2000 u*/ over W@ +
\ ( ymod caloff xscale -- )
	rot2
\ (xscale ymod caloff -- )
	dup 26 + W@ over 2+ W@ - rot d_2000 u*/ swap 2+ W@ +
\
	minY W@ max maxY W@ min swap
	minX W@ max maxX W@ min swap
;
\
\
: setXY
	6 x W@ y W@
	_scaledxy W@
	if
		_scalexy
	then
	5 swap
	in_debug W@ 8 and
	if
		c" x y -> sX sY: " .concstr x W@ .con y W@ .con 2 ST@ .con dup .con .concr
	then

	sm_setpos sm_setpos
;

: setMinMaxXY
	d_5000 dup dup dup minX W! minY W!  maxX W! maxY W!
	0 cal_table d_148 bounds

	do
		i W@    minX W@ over min minX W! maxX W@ max maxX W!
		i 2+ W@ minY W@ over min minY W! maxY W@ max maxY W!
	4 +loop
;

: cmd_proc
	4 state andnC!
	c" CMD_PROC" cds W!
	_red_led dup pinout pinhi
	_green_led dup pinout pinlo
	_yellow_led dup pinout pinlo
	_laser dup pinlo pinout
	0 in_topDown W!
	0 in_sideDown W!
	0 in_laserCount W!
	0 in_calibrate W!
	0 in_scan W!
	0 in_laserFlash W!
	_defClickIncrement _clickIncrement W!
	-1 _scaledxy W!
\
	d_5000 x W!
	d_5000 y W!
	d_102_000 d_170_000 5 sm_setminmax
	d_102_000 d_170_000 6 sm_setminmax
\
	setMinMaxXY
	d_5000 dup x W! y W! setXY
;
\
\
: topDown -1 in_topDown W! d_120 in_laserCount W! ;
: topUp 0 in_topDown W! ;
: sideDown -1 in_sideDown W! 0 in_laserCount W! ;
: sideUp 0 in_sideDown W! ;
\
: tick
	bootTicks L@ 1+ bootTicks L!
\
	in_laserCount W@ 0>
	in_laserFlash W@ dup 0= swap
	bootTicks L@ 1 and
\
	and or and _laser px
\
	in_laserCount W@ 1- 0 max in_laserCount W!
\
	bootTicks L@ cal_time <= _green_led px
\
	in_topDown W@ in_sideDown W@ and
	if
		bootTicks L@ cal_time <=
		if
			in_calibrate W@ 0= in_scan W@ 0= and
			if
				-1 in_calibrate W!
				lockdict c" calibrate" nfcog cogx  freedict
			then
		else
			in_calibrate W@ 0= in_scan W@ 0= and
			if
				-1 in_scan W!
				lockdict c" scan" nfcog cogx  freedict
			then
		then
	then
\
	in_calibrate W@ in_scan W@ or _yellow_led px

	in_debug W@ h10 and bootTicks L@ h1F and 0= and
	if
		_dbg
	then
;
\
: dbg in_debug W! ;
\
: topCW y W@ _clickIncrement W@ + d_10000 min y W! setXY ;
: topCCW y W@ _clickIncrement W@ - 0 max y W! setXY ;
\
: sideCW x W@ _clickIncrement W@ + d_10000 min x W! setXY ;
: sideCCW x W@ _clickIncrement W@ - 0 max x W! setXY ;
\
: scan
	4 state andnC!
	c" SCAN" cds W!
	c" Scan started~h0D~h0D" .concstr
	d_240 in_laserCount W!
\
	d_5000 0 do
		i 2* x W! i y W! setXY
		d_50 delms
	d_100 +loop
\
	d_5000 0 do
		d_10000 i 2* - x W! i d_5000 + y W! setXY
		d_50 delms
	d_100 +loop
\
	d_10000 0 do
		d_10000 i - y W! setXY
		d_50 delms
	d_100 +loop
\
	d_10000 0 do
		i x W! setXY
		d_50 delms
	d_100 +loop
\
	d_5000 0 do
		d_10000 i 2* - x W! i y W! setXY
		d_50 delms
	d_100 +loop
\
	d_5000 0 do
		i 2* x W! i d_5000 + y W! setXY
		d_50 delms
	d_100 +loop
\
	d_5000 0 do
		d_10000 i - y W! setXY
		d_50 delms
	d100 +loop
\
	d_5000 0 do
		d_10000 i - x W! setXY
		d_50 delms
	d100 +loop
\
	d_5000 dup x W! y W! setXY d500 delms
\
	0 in_laserCount W!
	0 in_scan W!
	4 state orC!
	0 cds W!
	c" Scan done~h0D~h0D" .concstr
;
\
: calibrate
	4 state andnC!
	c" CALIBRATE" cds W!
	-1 in_calibrate W!
	d_50 _clickIncrement W!
	0 _scaledxy  W!
	begin
		in_topDown W@ in_sideDown W@ and 0=
	until
\
	c" Calibration started~h0D~h0D" .concstr
\
	0 cal_table d_148 bounds
\
	do
		i W@ x W! i 2+ W@ y W! setXY
		c" Original point: " .concstr x W@ .con y W@ .con .concr
		d_240 in_laserCount W!
		d_16000 0
		do
			d_100 delms
			in_sideDown W@
			if
				in_topDown W@
				if
					drop -1
					begin
						in_sideDown W@ in_topDown W@ and 0=
					until
				else
					begin
						in_sideDown W@ 0=
					until
				then
\
				leave
			then
		loop
\
		dup if
			leave
		else
			x W@ i W!
			y W@ i 2+ W!
			c" Calibrated point: " .concstr x W@ .con y W@ .con .concr
		then
	4 +loop
\
	_defClickIncrement _clickIncrement W!
	-1 _scaledxy W!
\
	0= if
		d_5000 dup x W! y W! setXY
		d_40 in_laserCount W!
		-1 in_laserFlash W!
		c" Calibration write?~h0D~h0D" .concstr
\
		d_1000 0
		do
			d_100 delms
			in_topDown W@ in_sideDown W@ and
			if
				c" Calibration write started~h0D~h0D" .concstr
				0 in_laserFlash W!
				lockdict cogid nfcog iolink freedict
				." fsdrop~h0D~h0D"
				." fswrite cal.f~h0D~h0D"
				." lockdict wvariable cal_table -2 allot" cr
				0 cal_table d_148 bounds
				do
					." d_" i W@ . ." w, "
					1+ dup 12 >=
					if
						drop 0 cr
					then
				2 +loop
				drop
				." freedict~h0D~h0D...~h0D~h0D"
				cogid iounlink
				c" Calibration write done~h0D~h0D" .concstr
\
				d_50 in_laserCount W!
				d_5000 delms
\
				leave
			then
		loop
	then
\
	0 in_laserFlash W!
	0 in_laserCount W!
\
	0 in_calibrate W!
	4 state orC!
	0 cds W!
	c" Calibration done~h0D~h0D" .concstr
;
\
\
c" 848-" prop W@ ccopy
: onreset6
	fkey? and fkey? and or h1B <>
	if
		sm_start_servos
		c" in_monitor" 0 cogx
		c" cmd_proc" 2 cogx
		0 2 iolink
\
		$S_con iodis $S_con cogreset 100 delms
		c" hcRx hcTx 57600 serial" $S_con cogx 100 delms
		cogid >con
	then
	sbv .concstr
	c" onreset6" (forget)
;
\
\

...

\
fswrite cal.f

lockdict wvariable cal_table -2 allot
d_2950 w, d_900 w, d_4150 w, d_850 w, d_5450 w, d_750 w, d_6850 w, d_750 w, d_8250 w, d_950 w, d_9400 w, d_1150 w,
d_2800 w, d_2150 w, d_4050 w, d_2100 w, d_5350 w, d_2100 w, d_6800 w, d_2150 w, d_8200 w, d_2250 w, d_9300 w, d_2350 w,
d_2700 w, d_3600 w, d_3850 w, d_3550 w, d_5250 w, d_3550 w, d_6800 w, d_3600 w, d_8100 w, d_3650 w, d_9300 w, d_3750 w,
d_2650 w, d_4950 w, d_3800 w, d_5050 w, d_5300 w, d_5100 w, d_6700 w, d_5050 w, d_8000 w, d_5150 w, d_9250 w, d_5100 w,
d_2500 w, d_6350 w, d_3700 w, d_6500 w, d_5100 w, d_6550 w, d_6600 w, d_6600 w, d_7950 w, d_6600 w, d_9050 w, d_6550 w,
d_2450 w, d_7650 w, d_3600 w, d_7850 w, d_6250 w, d_8050 w, d_6550 w, d_7900 w, d_7900 w, d_7900 w, d_9100 w, d_7800 w,
d_6042 w, d_4367 w, freedict

...

ext: f


text:
! { dg-do run }
! { dg-add-options ieee }
! { dg-skip-if "NaN not supported" { spu-*-* } { "*" } { "" } }
!
! PR fortran/34209
!
! Test run-time implementation of NEAREST
!
program test
  implicit none
  real(4), volatile :: r4
  real(8), volatile :: r8

! Single precision with single-precision sign

  r4 = 0.0_4
  ! 0+ > 0
  if (nearest(r4, 1.0) &
      <= r4) &
    call abort()
  ! 0++ > 0+
  if (nearest(nearest(r4, 1.0), 1.0) &
      <= nearest(r4, 1.0)) &
    call abort()
  ! 0+++ > 0++
  if (nearest(nearest(nearest(r4, 1.0), 1.0), 1.0) &
      <= nearest(nearest(r4, 1.0), 1.0)) &
    call abort()
  ! 0+- = 0
  if (nearest(nearest(r4, 1.0), -1.0) &
      /= r4) &
    call abort()
  ! 0++- = 0+
  if (nearest(nearest(nearest(r4, 1.0), 1.0), -1.0) &
      /= nearest(r4, 1.0)) &
    call abort()
  ! 0++-- = 0
  if (nearest(nearest(nearest(nearest(r4, 1.0), 1.0), -1.0), -1.0) &
      /= r4) &
    call abort()

  ! 0- < 0
  if (nearest(r4, -1.0) &
      >= r4) &
    call abort()
  ! 0-- < 0+
  if (nearest(nearest(r4, -1.0), -1.0) &
      >= nearest(r4, -1.0)) &
    call abort()
  ! 0--- < 0--
  if (nearest(nearest(nearest(r4, -1.0), -1.0), -1.0) &
      >= nearest(nearest(r4, -1.0), -1.0)) &
    call abort()
  ! 0-+ = 0
  if (nearest(nearest(r4, -1.0), 1.0) &
      /= r4) &
    call abort()
  ! 0--+ = 0-
  if (nearest(nearest(nearest(r4, -1.0), -1.0), 1.0) &
      /= nearest(r4, -1.0)) &
    call abort()
  ! 0--++ = 0
  if (nearest(nearest(nearest(nearest(r4, -1.0), -1.0), 1.0), 1.0) &
      /= r4) &
    call abort()

  r4 = 42.0_4
  ! 42++ > 42+
  if (nearest(nearest(r4, 1.0), 1.0) &
      <= nearest(r4, 1.0)) &
    call abort()
  ! 42-- < 42-
  if (nearest(nearest(r4, -1.0), -1.0) &
      >= nearest(r4, -1.0)) &
    call abort()
  ! 42-+ = 42
  if (nearest(nearest(r4, -1.0), 1.0) &
      /= r4) &
    call abort()
  ! 42+- = 42
  if (nearest(nearest(r4, 1.0), -1.0) &
      /= r4) &
    call abort()

  r4 = 0.0
  ! INF+ = INF
  if (nearest(1.0/r4, 1.0) /= 1.0/r4) call abort()
  ! -INF- = -INF
  if (nearest(-1.0/r4, -1.0) /= -1.0/r4) call abort()
  ! NAN- = NAN
  if (.not.isnan(nearest(0.0/r4,  1.0))) call abort()
  ! NAN+ = NAN
  if (.not.isnan(nearest(0.0/r4, -1.0))) call abort()

! Double precision with single-precision sign

  r8 = 0.0_8
  ! 0+ > 0
  if (nearest(r8, 1.0) &
      <= r8) &
    call abort()
  ! 0++ > 0+
  if (nearest(nearest(r8, 1.0), 1.0) &
      <= nearest(r8, 1.0)) &
    call abort()
  ! 0+++ > 0++
  if (nearest(nearest(nearest(r8, 1.0), 1.0), 1.0) &
      <= nearest(nearest(r8, 1.0), 1.0)) &
    call abort()
  ! 0+- = 0
  if (nearest(nearest(r8, 1.0), -1.0) &
      /= r8) &
    call abort()
  ! 0++- = 0+
  if (nearest(nearest(nearest(r8, 1.0), 1.0), -1.0) &
      /= nearest(r8, 1.0)) &
    call abort()
  ! 0++-- = 0
  if (nearest(nearest(nearest(nearest(r8, 1.0), 1.0), -1.0), -1.0) &
      /= r8) &
    call abort()

  ! 0- < 0
  if (nearest(r8, -1.0) &
      >= r8) &
    call abort()
  ! 0-- < 0+
  if (nearest(nearest(r8, -1.0), -1.0) &
      >= nearest(r8, -1.0)) &
    call abort()
  ! 0--- < 0--
  if (nearest(nearest(nearest(r8, -1.0), -1.0), -1.0) &
      >= nearest(nearest(r8, -1.0), -1.0)) &
    call abort()
  ! 0-+ = 0
  if (nearest(nearest(r8, -1.0), 1.0) &
      /= r8) &
    call abort()
  ! 0--+ = 0-
  if (nearest(nearest(nearest(r8, -1.0), -1.0), 1.0) &
      /= nearest(r8, -1.0)) &
    call abort()
  ! 0--++ = 0
  if (nearest(nearest(nearest(nearest(r8, -1.0), -1.0), 1.0), 1.0) &
      /= r8) &
    call abort()

  r8 = 42.0_8
  ! 42++ > 42+
  if (nearest(nearest(r8, 1.0), 1.0) &
      <= nearest(r8, 1.0)) &
    call abort()
  ! 42-- < 42-
  if (nearest(nearest(r8, -1.0), -1.0) &
      >= nearest(r8, -1.0)) &
    call abort()
  ! 42-+ = 42
  if (nearest(nearest(r8, -1.0), 1.0) &
      /= r8) &
    call abort()
  ! 42+- = 42
  if (nearest(nearest(r8, 1.0), -1.0) &
      /= r8) &
    call abort()

  r4 = 0.0
  ! INF+ = INF
  if (nearest(1.0/r4, 1.0) /= 1.0/r4) call abort()
  ! -INF- = -INF
  if (nearest(-1.0/r4, -1.0) /= -1.0/r4) call abort()
  ! NAN- = NAN
  if (.not.isnan(nearest(0.0/r4,  1.0))) call abort()
  ! NAN+ = NAN
  if (.not.isnan(nearest(0.0/r4, -1.0))) call abort()

! Single precision with double-precision sign

  r4 = 0.0_4
  ! 0+ > 0
  if (nearest(r4, 1.0d0) &
      <= r4) &
    call abort()
  ! 0++ > 0+
  if (nearest(nearest(r4, 1.0d0), 1.0d0) &
      <= nearest(r4, 1.0d0)) &
    call abort()
  ! 0+++ > 0++
  if (nearest(nearest(nearest(r4, 1.0d0), 1.0d0), 1.0d0) &
      <= nearest(nearest(r4, 1.0d0), 1.0d0)) &
    call abort()
  ! 0+- = 0
  if (nearest(nearest(r4, 1.0d0), -1.0d0) &
      /= r4) &
    call abort()
  ! 0++- = 0+
  if (nearest(nearest(nearest(r4, 1.0d0), 1.0d0), -1.0d0) &
      /= nearest(r4, 1.0d0)) &
    call abort()
  ! 0++-- = 0
  if (nearest(nearest(nearest(nearest(r4, 1.0d0), 1.0d0), -1.0d0), -1.0d0) &
      /= r4) &
    call abort()

  ! 0- < 0
  if (nearest(r4, -1.0d0) &
      >= r4) &
    call abort()
  ! 0-- < 0+
  if (nearest(nearest(r4, -1.0d0), -1.0d0) &
      >= nearest(r4, -1.0d0)) &
    call abort()
  ! 0--- < 0--
  if (nearest(nearest(nearest(r4, -1.0d0), -1.0d0), -1.0d0) &
      >= nearest(nearest(r4, -1.0d0), -1.0d0)) &
    call abort()
  ! 0-+ = 0
  if (nearest(nearest(r4, -1.0d0), 1.0d0) &
      /= r4) &
    call abort()
  ! 0--+ = 0-
  if (nearest(nearest(nearest(r4, -1.0d0), -1.0d0), 1.0d0) &
      /= nearest(r4, -1.0d0)) &
    call abort()
  ! 0--++ = 0
  if (nearest(nearest(nearest(nearest(r4, -1.0d0), -1.0d0), 1.0d0), 1.0d0) &
      /= r4) &
    call abort()

  r4 = 42.0_4
  ! 42++ > 42+
  if (nearest(nearest(r4, 1.0d0), 1.0d0) &
      <= nearest(r4, 1.0d0)) &
    call abort()
  ! 42-- < 42-
  if (nearest(nearest(r4, -1.0d0), -1.0d0) &
      >= nearest(r4, -1.0d0)) &
    call abort()
  ! 42-+ = 42
  if (nearest(nearest(r4, -1.0d0), 1.0d0) &
      /= r4) &
    call abort()
  ! 42+- = 42
  if (nearest(nearest(r4, 1.0d0), -1.0d0) &
      /= r4) &
    call abort()

  r4 = 0.0
  ! INF+ = INF
  if (nearest(1.0d0/r4, 1.0d0) /= 1.0d0/r4) call abort()
  ! -INF- = -INF
  if (nearest(-1.0d0/r4, -1.0d0) /= -1.0d0/r4) call abort()
  ! NAN- = NAN
  if (.not.isnan(nearest(0.0/r4,  1.0d0))) call abort()
  ! NAN+ = NAN
  if (.not.isnan(nearest(0.0/r4, -1.0d0))) call abort()

! Double precision with double-precision sign

  r8 = 0.0_8
  ! 0+ > 0
  if (nearest(r8, 1.0d0) &
      <= r8) &
    call abort()
  ! 0++ > 0+
  if (nearest(nearest(r8, 1.0d0), 1.0d0) &
      <= nearest(r8, 1.0d0)) &
    call abort()
  ! 0+++ > 0++
  if (nearest(nearest(nearest(r8, 1.0d0), 1.0d0), 1.0d0) &
      <= nearest(nearest(r8, 1.0d0), 1.0d0)) &
    call abort()
  ! 0+- = 0
  if (nearest(nearest(r8, 1.0d0), -1.0d0) &
      /= r8) &
    call abort()
  ! 0++- = 0+
  if (nearest(nearest(nearest(r8, 1.0d0), 1.0d0), -1.0d0) &
      /= nearest(r8, 1.0d0)) &
    call abort()
  ! 0++-- = 0
  if (nearest(nearest(nearest(nearest(r8, 1.0d0), 1.0d0), -1.0d0), -1.0d0) &
      /= r8) &
    call abort()

  ! 0- < 0
  if (nearest(r8, -1.0d0) &
      >= r8) &
    call abort()
  ! 0-- < 0+
  if (nearest(nearest(r8, -1.0d0), -1.0d0) &
      >= nearest(r8, -1.0d0)) &
    call abort()
  ! 0--- < 0--
  if (nearest(nearest(nearest(r8, -1.0d0), -1.0d0), -1.0d0) &
      >= nearest(nearest(r8, -1.0d0), -1.0d0)) &
    call abort()
  ! 0-+ = 0
  if (nearest(nearest(r8, -1.0d0), 1.0d0) &
      /= r8) &
    call abort()
  ! 0--+ = 0-
  if (nearest(nearest(nearest(r8, -1.0d0), -1.0d0), 1.0d0) &
      /= nearest(r8, -1.0d0)) &
    call abort()
  ! 0--++ = 0
  if (nearest(nearest(nearest(nearest(r8, -1.0d0), -1.0d0), 1.0d0), 1.0d0) &
      /= r8) &
    call abort()

  r8 = 42.0_8
  ! 42++ > 42+
  if (nearest(nearest(r8, 1.0d0), 1.0d0) &
      <= nearest(r8, 1.0d0)) &
    call abort()
  ! 42-- < 42-
  if (nearest(nearest(r8, -1.0d0), -1.0d0) &
      >= nearest(r8, -1.0d0)) &
    call abort()
  ! 42-+ = 42
  if (nearest(nearest(r8, -1.0d0), 1.0d0) &
      /= r8) &
    call abort()
  ! 42+- = 42
  if (nearest(nearest(r8, 1.0d0), -1.0d0) &
      /= r8) &
    call abort()

  r4 = 0.0
  ! INF+ = INF
  if (nearest(1.0d0/r4, 1.0d0) /= 1.0d0/r4) call abort()
  ! -INF- = -INF
  if (nearest(-1.0d0/r4, -1.0d0) /= -1.0d0/r4) call abort()
  ! NAN- = NAN
  if (.not.isnan(nearest(0.0/r4,  1.0d0))) call abort()
  ! NAN+ = NAN
  if (.not.isnan(nearest(0.0/r4, -1.0d0))) call abort()

end program test

max_stars_repo_name: maurizioabba/rose

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


text:
      SUBROUTINE strip_comments(input_file,output_file)
      USE safe_open_mod
      IMPLICIT NONE
!
!     strips comment lines (starting with '!') from input_file
!     renames clean file input_file // '.stripped'
!
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      CHARACTER(LEN=*) :: input_file
      CHARACTER(LEN=*) :: output_file
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER, PARAMETER :: unit_strip = 20
      INTEGER :: istat, iustrip, iunew, itry
      CHARACTER(LEN=1024) :: line
      LOGICAL :: lex
      INTEGER, EXTERNAL :: getcwd

C-----------------------------------------------
      itry = 1
 99   INQUIRE (file=input_file, exist=lex, iostat=istat)
      IF (.not.lex .and. itry .lt. 200) THEN
         itry = itry + 1
         GOTO 99
      END IF
      IF (istat.ne.0 .or. .not. lex) THEN
!         istat = getcwd(line)
         PRINT *,TRIM(input_file),' does not exist in directory: ',
     1           TRIM(line)
         CALL FLUSH(6)
         STOP
      END IF

      iustrip = unit_strip
      CALL safe_open(iustrip, istat, input_file, 'old','formatted')
      IF (istat .ne. 0) THEN
        line = 'error opening ' // TRIM(input_file) //
     1         ' in STRIP_COMMENTS'
        PRINT *, line
        PRINT *,'istat = ', istat
        STOP
      END IF
      iunew = iustrip + 1
      CALL safe_open(iunew, istat, TRIM(output_file),
     1     'replace', 'formatted')
!      CALL safe_open(iunew, istat, TRIM(input_file) // '.stripped',
!     1     'replace', 'formatted')
      IF (istat .ne. 0) THEN
        line = 'error opening ' // TRIM(output_file) //
     1      ' in STRIP_COMMENTS'
!        line = 'error opening ' // TRIM(input_file) //
!     1      '.stripped  in STRIP_COMMENTS'
        PRINT *, line
        PRINT *,'istat = ', istat
        STOP
      END IF
      DO
         READ(iustrip, '(a)', END=100) line
         line = ADJUSTL(line)
         IF (line(1:1) == '!') CYCLE
         WRITE(iunew, '(a)') TRIM(line)
      END DO

      CALL FLUSH(iunew)

 100  CONTINUE

      CLOSE(iustrip)
      CLOSE(iunew)

      END SUBROUTINE strip_comments

max_stars_repo_name: joseluisvelasco/STELLOPT

max_issues_repo_name: joseluisvelasco/STELLOPT

max_forks_repo_name: joseluisvelasco/STELLOPT


text:
Conference and Event Services (formerly called Campus Events & Visitor Services)

Planning an event? Look no further than Conference and Event Services at UC Davis. Our event planning team specializes in corporate, athletic, and special events, and offers logistical planning, onsite coordination, and all the amenities of the University. The distinct combination of charm and prestige our campus has made hundreds of guests choose UC Davis as their venue to say I do, celebrate retirement, and hold family reunions. Conference and Event Services has over 15 years of event planning experience as Davis premier event coordinator.

Hotel and Conference Center     198533

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
prime_do_while/reset.cpp
prime_do_while/display.cpp
prime_do_while/prime_numgen.cpp
prime_do_while/main.cpp


! Test specification expressions

module m
  type :: t(n)
    integer, len :: n = 1
    character(len=n) :: c
  end type
  interface
    integer function foo()
    end function
    pure real function realfunc(x)
      real, intent(in) :: x
    end function
    pure integer function hasProcArg(p)
      import realfunc
      procedure(realfunc) :: p
    end function
  end interface
  integer :: coarray[*]
 contains
  pure integer function modulefunc1(n)
    integer, value :: n
    modulefunc1 = n
  end function
  subroutine test(out, optional)
    !ERROR: Invalid specification expression: reference to impure function 'foo'
    type(t(foo())) :: x1
    integer :: local
    !ERROR: Invalid specification expression: reference to local entity 'local'
    type(t(local)) :: x2
    !ERROR: The internal function 'internal' cannot be referenced in a specification expression
    type(t(internal(0))) :: x3
    integer, intent(out) :: out
    !ERROR: Invalid specification expression: reference to INTENT(OUT) dummy argument 'out'
    type(t(out)) :: x4
    integer, intent(in), optional :: optional
    !ERROR: Invalid specification expression: reference to OPTIONAL dummy argument 'optional'
    type(t(optional)) :: x5
    !ERROR: Invalid specification expression: dummy procedure argument
    type(t(hasProcArg(realfunc))) :: x6
    !ERROR: Invalid specification expression: coindexed reference
    type(t(coarray[1])) :: x7
    type(t(kind(foo()))) :: x101 ! ok
    type(t(modulefunc1(0))) :: x102 ! ok
    !ERROR: The module function 'modulefunc2' must have been previously defined when referenced in a specification expression
    type(t(modulefunc2(0))) :: x103 ! ok
   contains
    pure integer function internal(n)
      integer, value :: n
      internal = n
    end function
  end subroutine
  pure integer function modulefunc2(n)
    integer, value :: n
    modulefunc2 = n
  end function
end module


text:
Ross Dress for Less opened in Woodland 20051001. Ross is a discount clothing store similar to Marshalls and TJ Maxx. They sell overstock or slightly irregular brand name clothes at a heavy discount compared to full department store price. Next door to The Tractor Supply Co.

A Ross was originally going to be built in Davis years ago but was shut out by the Davis City Council for not being upscale enough, and the oft repeated fear that it would threaten downtown business.  There was a bit of an uproar about it since at the time there were few, if any, affordable clothing stores in town.  It was said at the time that you couldnt buy a pair of underwear in town.

Here is the inside story from a former Ross employee: The brand name clothes are cheap because the clothing that Ross sells is either overstock or slightly irregular. Ross gets cheap clothes that other retailers were not able to sell or items that are defective. Examples of the stores where the clothes have come from are Macys and Nordstrom. Many of the clothes on sale have been worn, tried on, and returned to the store over and over due to Ross lax return policy. The clothes are not cleaned prior to being put back on the racks. Please shower as soon as possible after trying on clothes at a Ross Store.

Ross does not actually have Sales in a traditional sense. In essence, they have an everyday low price policy. The prices are not lowered one week and then raised like at Macys. Items are at a set price and then go down the longer they sit on the clothing racks. So the price you see on the item is the highest price and it may ring up as less at the register.

Ross also offers housewares, dishes, pots, pans, utensils, chocolate, candies, fine jewelry, expensive perfumes at a heavy discount, designer sunglasses, books, lamps, chairs, toys, sports equipment, shoes, small electric devices, luggage, purses, and a limited selection of music and movies.

To get to Ross from Davis, take Yolobus 42A to E. Main at Matmor.


text:
      subroutine find_closest_curve(flag, par1v, par2v, par3v,
     2     i, j, k, points, weight)
!***********************************************************************
!  Copyright, 2004,  The  Regents  of the  University of California.
!  This program was prepared by the Regents of the University of
!  California at Los Alamos National Laboratory (the University) under
!  contract No. W-7405-ENG-36 with the U.S. Department of Energy (DOE).
!  All rights in the program are reserved by the DOE and the University.
!  Permission is granted to the public to copy and use this software
!  without charge, provided that this Notice and any statement of
!  authorship are reproduced on all copies. Neither the U.S. Government
!  nor the University makes any warranty, express or implied, or
!  assumes any liability or responsibility for the use of this software.
!**********************************************************************
!D1
!D1 PURPOSE
!D1
!D1 To determine the curve with the closest parameter values to the input
!D1
!**********************************************************************
!D2
!D2 REVISION HISTORY
!D2
!D2 FEHM Version 2.20 [10086-STN-2.20-00]
!D2
!D2 Initial implementation: 15-DEC-02, Programmer: Bruce Robinson
!D2
!D2 $Log:   /pvcs.config/fehm90/src/find_closest_curve.f_a  $
!D2
!D2    Rev 2.5   06 Jan 2004 10:43:00   pvcs
!D2 FEHM Version 2.21, STN 10086-2.21-00, Qualified October 2003
!D2
!D2    Rev 2.4   29 Jan 2003 09:21:48   pvcs
!D2 FEHM Version 2.20, STN 10086-2.20-00
!D2
!**********************************************************************
!D3
!D3 REQUIREMENTS TRACEABILITY
!D3
!D3 2.3.5 Cell-based particle-tracking module
!D3 2.3.6 Streamline particle-tracking module
!D3
!**********************************************************************
!D4
!D4 SPECIAL COMMENTS AND REFERENCES
!D4
!D4 Requirements from SDN: 10086-RD-2.20-00
!D4   SOFTWARE REQUIREMENTS DOCUMENT (RD) for the
!D4   FEHM Application Version 2.20
!D4
!**********************************************************************
! Determine the transfer function curve with parameters closest to
! the input values

      use compfrac
      use comai, only : ierr
      implicit none

      real*8 par1v, par2v, par3v
      real*8 pl1, pl2, pl3
      integer i, j, k
      integer i1
      real*8, allocatable, dimension(:) :: distance
      real*8, allocatable, dimension(:) :: p1
      real*8, allocatable, dimension(:) :: p2
      real*8, allocatable, dimension(:) :: p3
      integer, allocatable, dimension(:) :: iorder
      integer flag
      integer nparamp1, ifind, icol
      integer points(4)
      real weight(4), a(4,4), v(4,4), w(4), wttot
      integer ioindex(4), return_flag
      return_flag = 0

      allocate(distance(nump1), p1(nump1),p2(nump1),p3(nump1))
      allocate(iorder(nump1))
      distance = 0.
      iorder = 0
      p1 = 0.
      p2 = 0.
      p3 = 0.
      nparamp1=numparams+1

c     Load the parameter arrays and convert the appropriate
c     values to log values where specified
      if(log_flag(1).eq.1) then
         pl1 = dlog(par1v)/alog(normal_param(1))
         do i1 = 1, nump1
            p1(i1) = dlog(param1(i1))/alog(normal_param(1))
         end do
      else
         pl1 = par1v/normal_param(1)
         do i1 = 1, nump1
            p1(i1) = param1(i1)/normal_param(1)
         end do
      end if
      if(log_flag(2).eq.1) then
         pl2 = dlog(par2v)/alog(normal_param(2))
         do i1 = 1, nump1
            p2(i1) = dlog(param2(i1))/alog(normal_param(2))
         end do
      else
         pl2 = par2v/normal_param(2)
         do i1 = 1, nump1
            p2(i1) = param2(i1)/normal_param(2)
         end do
      end if
      if(numparams.gt.2) then
         if(log_flag(3).eq.1) then
            pl3 = dlog(par3v)/alog(normal_param(3))
            do i1 = 1, nump1
               p3(i1) = dlog(param3(i1))/alog(normal_param(3))
            end do
         else
            pl3 = par3v/normal_param(3)
            do i1 = 1, nump1
               p3(i1) = param3(i1)/normal_param(3)
            end do
         end if
      else
         pl3 = par3v
         do i1 = 1, nump1
            p3(i1) = pl3
         end do
      end if

c     Compute distance from the parameter to each curve

      do i1 = 1, nump1
         distance(i1) = sqrt((p1(i1)-pl1)**2+
     2        (p2(i1)-pl2)**2+
     3        (p3(i1)-pl3)**2)
      end do
c     Sort from lowest value to highest
      call sort_values(nump1, nump1, iorder, distance)

      if(flag.eq.0) then

c     Index of curve that has smallest deviation from the
c     input curve. j and k are 1 because of the data structure
c     i.e. curves go from 1 to nump1 in i, and j and k are 1 long

         i = iorder(1)
         j = 1
         k = 1

      else

c     Determine if any of the parameters are identical over all
c     points and replace with points somewhat farther away.
c     Always use closest two points (unless they are identical)
c     Start replacing 3rd and 4th points as necessary to get
c     some differences in each parameter (otherwise SVD
c     fails)

         call  filter_points

c     Set parameters to do nearest neighbor if an appropriate
c     set of points could not be found for svdcmp

         if(return_flag.ne.0) then
            if(numparams.eq.2) then
               weight(1) = 1.
               weight(2) = 0.
               weight(3) = 0.
               points(1) = iorder(1)
               points(2) = iorder(2)
               points(3) = iorder(3)
            elseif(numparams.gt.2) then
               weight(1) = 1.
               weight(2) = 0.
               weight(3) = 0.
               weight(4) = 0.
               points(1) = iorder(1)
               points(2) = iorder(2)
               points(3) = iorder(3)
               points(4) = iorder(4)
            end if
            return
         end if

c     Least squares fit to numparams+1 nearest curves
c     using singular value decomposition to find weights
c     to be used in interpolation
         a = 0.
         a(1,1) = pl1-p1(iorder(ioindex(1)))
         a(1,2) = pl1-p1(iorder(ioindex(2)))
         a(1,3) = pl1-p1(iorder(ioindex(3)))
         a(2,1) = pl2-p2(iorder(ioindex(1)))
         a(2,2) = pl2-p2(iorder(ioindex(2)))
         a(2,3) = pl2-p2(iorder(ioindex(3)))
         if(numparams.eq.3) then
            a(1,4) = pl1-p1(iorder(ioindex(4)))
            a(2,4) = pl2-p2(iorder(ioindex(4)))
            a(3,1) = pl3-p3(iorder(ioindex(1)))
            a(3,2) = pl3-p3(iorder(ioindex(2)))
            a(3,3) = pl3-p3(iorder(ioindex(3)))
            a(3,4) = pl3-p3(iorder(ioindex(4)))
         end if

c     SVD algorithm finds the weights

         call svdcmp_new(ierr,a,nparamp1,nparamp1,4,4,w,v)

         ifind = 0
         find_soln: do i = 1, nparamp1
            if(w(i).lt.1.e-6) then
               ifind = i
               exit find_soln
            end if
         end do find_soln

         if(ifind.eq.0) then
            if(numparams.eq.2) then
               weight(1) = 1.
               weight(2) = 0.
               weight(3) = 0.
               points(1) = iorder(1)
               points(2) = iorder(2)
               points(3) = iorder(3)
               write(ierr,*) 'Error in particle tracking interp.'
               write(ierr,*) 'Stopping'
            stop
            elseif(numparams.gt.2) then
               weight(1) = 1.
               weight(2) = 0.
               weight(3) = 0.
               weight(4) = 0.
               points(1) = iorder(1)
               points(2) = iorder(2)
               points(3) = iorder(3)
               points(4) = iorder(4)
            end if
            return
         end if

         weight = 0.
         wttot = 0.
         do icol = 1, nparamp1
            weight(icol) = v(icol,ifind)
            wttot = wttot + weight(icol)
         end do

c     Check and revert to nearest neighbor if poorly behaved

         do icol = 1, nparamp1
            if(abs(weight(icol)).gt.1e-10) then
c     weight_factor is a term that prevents the svd results to be used
c     when the weights themselves are much larger than the sum
c     of the weights. This is a situation in which numerical errors
c     can become a problem in the interpolation scheme. The code
c     reverts to the nearest neighbor if this happens.

               if(abs(wttot/weight(icol)).lt.weight_factor) then
                  wttot = 1.
                  weight(1) = 1.
                  weight(2) = 0.
                  weight(3) = 0.
                  if(numparams.gt.2) then
                     weight(4) = 0.
                  end if
               end if
            end if
         end do
         weight = weight/wttot
         points(1) = iorder(ioindex(1))
         points(2) = iorder(ioindex(2))
         points(3) = iorder(ioindex(3))
         if(numparams.gt.2) then
            points(4) = iorder(ioindex(4))
         end if
      end if

      deallocate(distance, iorder, p1, p2, p3)

      return
      contains

      subroutine filter_points
      implicit none
      real min_difference
      parameter (min_difference = 1.e-3)
      integer ipflag(3), ipts, jpar
      real pmatrix(3,4), denom, comp_value
      real max_p(3)
      logical done

      ioindex(1) = 1
      ioindex(2) = 2
      ioindex(3) = 3
      ioindex(4) = 4

      done = .false.
      do while(.not.done)
         ipflag = 0
c     Load matrix with inital values for checking the parameter values
         max_p = 0.
         do ipts = 1, nparamp1
            pmatrix(1,ipts) = p1(iorder(ioindex(ipts)))
            pmatrix(2,ipts) = p2(iorder(ioindex(ipts)))
            if(numparams.gt.2) then
               pmatrix(3,ipts) = p3(iorder(ioindex(ipts)))
            end if
         end do

         do jpar = 1, numparams

            inner_loop: do ipts = 2, nparamp1
               denom = max(abs(pmatrix(jpar,ipts)),
     2           abs(pmatrix(jpar,ipts-1)))
               if(denom.ne.0) then
                  comp_value = abs(pmatrix(jpar,ipts)-
     2                 pmatrix(jpar,ipts-1))/denom
                  if(comp_value.gt.min_difference) then
                     ipflag(jpar) = 1
                     exit inner_loop
                  else
                  end if
               else
               end if
            end do inner_loop

         end do

         done = .true.
         do jpar = 1, numparams
            if(ipflag(jpar).eq.0) done = .false.
         end do
         if(done) then
         end if
         if(.not.done) then
            if(numparams.eq.2) then
               if(ioindex(3).lt.nump1) then
                  ioindex(3) = ioindex(3) + 1
               else
                  return_flag = 1
                  done = .true.
               end if
            elseif(numparams.gt.2) then
               if(ioindex(4).lt.nump1) then
                  ioindex(3) = ioindex(3) + 1
                  ioindex(4) = ioindex(4) + 1
               else
                  return_flag = 1
                  done = .true.
               end if
            end if
         end if
      end do

      end subroutine filter_points
      end subroutine find_closest_curve

max_stars_repo_name: satkarra/FEHM

max_issues_repo_name: satkarra/FEHM

max_forks_repo_name: satkarra/FEHM


text:
** Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
** See https://llvm.org/LICENSE.txt for license information.
** SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception

*   Register allocation bug when passing arguments by value to a
*   STDCALL *   routine on win32.  Even though STDCALL is only available on
*   win32, this test can still be used on other systems; there will not be
*   any value passing, but that's ok.

        integer function ifoo( hInstance, nCmdShow)
!DEC$ ATTRIBUTES STDCALL :: ifoo
        integer hInstance
        integer nCmdShow
        ii = hinstance
        jj = ncmdshow
        call bar(ii,jj)
        ifoo = 0
        end
        subroutine bar(ii,jj)
	common/ires/ires(2)
	ires(1) = ii
	ires(2) = jj
        end
	common/ires/ires(2)
	integer iexp(2)
	data iexp/1,2/
        integer ifoo
        external ifoo
!DEC$ ATTRIBUTES STDCALL :: ifoo
        kk = ifoo(1,2)
	call check(ires, iexp, 2)
        end

max_stars_repo_name: abrahamtovarmob/flang

max_issues_repo_name: abrahamtovarmob/flang

max_forks_repo_name: abrahamtovarmob/flang


text:
      MODULE efsav_I
      INTERFACE
!...Generated by Pacific-Sierra Research 77to90  4.4G  10:47:10  03/09/06
      SUBROUTINE efsav (TT0, HESS, FUNCT, GRAD, XPARAM, PMAT, IL, JL, BMAT, IPOW)
      USE vast_kind_param,ONLY: DOUBLE
      doubleprecision DUMY
      real(DOUBLE), INTENT(INOUT) :: TT0
      real(DOUBLE), DIMENSION(*), INTENT(IN) :: HESS
      real(DOUBLE), INTENT(IN) :: FUNCT
      real(DOUBLE), DIMENSION(*), INTENT(IN) :: GRAD
      real(DOUBLE), DIMENSION(*), INTENT(IN) :: XPARAM
      real(DOUBLE), DIMENSION(*), INTENT(IN) :: PMAT
      integer, INTENT(IN) :: IL
      integer, INTENT(IN) :: JL
      real(DOUBLE), DIMENSION(*), INTENT(IN) :: BMAT
      integer, DIMENSION(9), INTENT(INOUT) :: IPOW
      END SUBROUTINE
      END INTERFACE
      END MODULE

max_stars_repo_name: openmopac/MOPAC-archive

max_issues_repo_name: openmopac/MOPAC-archive

max_forks_repo_name: openmopac/MOPAC-archive


text:
logical function isanynan(arr)
    use params, only : dp

    implicit none

    real(dp), dimension(:,:) :: arr
    integer, dimension(2) :: shapes
    integer :: i, j

    shapes = shape(arr)
    isanynan = .false.
    do i = 1, shapes(1)
        do j = 1, shapes(2)
            if (isnan(arr(i,j))) then
                isanynan = .true.
                return
            endif
        enddo
    enddo
endfunction

logical function b_any(arr)
    implicit none
    logical :: arr(:,:)
    b_any = any(arr)
end function
    ! integer :: r

    ! r = rank(arr)
    ! isanynan = isanynan_helper(arr, r)
! end function

! logical isanynan_helper(arr, rank)
    ! real(REAL64), intent(in), dimension(*) :: arr
    ! integer, intent(in) :: rank
    ! integer, dimension(rank) :: sizes

    ! sizes = shape(arr)
    ! do i = 1,

max_stars_repo_name: SpakowitzLab/BasicWLC

max_issues_repo_name: riscalab/wlcsim

max_forks_repo_name: riscalab/wlcsim


text:
Program Ejemplo18

!
! Uso de Format para escribir
!

IMPLICIT NONE

Character(4) :: Char

Integer ii, jj

Real :: xx, yy, uu

Real(8) :: zz, ww

Char = 'Cuat'

uu = 1

ii = 12345678; jj = 34

Print *, 'escrito generico :',ii

xx = 2.345678901234; yy = 0.012345

zz = 1.234567890123456789D0

ww = 1.0D0

 zz = zz + 1

!  zz = 1.234567890123456789D0

! zz = 1.234567890123456789

Print 10, ii

10 Format(I9)

100 Format(// 11X, I4, 2X, I5, 4X, F22.14, 2X, G22.14E4, 2X, E22.14E4)

Print 100, ii, jj, yy, yy, yy

Print 100, jj, jj, zz, zz, zz
Print 100, jj

Print 200, zz, Char

200 Format(// 5X, E22.14, 5X, A5)

Print 200, zz

Print 300, zz

300 Format(5X, E35.27)

Print *, 'ww :', ww
Print 300, ww
Print 200, ww, Char

Print*,'********************************'

400 Format(// 5X,'****', 3(3X,I10))

print 400, ii, jj, ii

End Program Ejemplo18

max_stars_repo_name: EdgardoBonzi/Fortran-Examples

max_issues_repo_name: EdgardoBonzi/Fortran-Examples

max_forks_repo_name: EdgardoBonzi/Fortran-Examples


text:
!! R1001 format-stmt
!     is  FORMAT format-specification
!
!! R1002 format-specification
!     is  ( [ format-item-list ] )
!     or  ( [ format-items, ] unlimited-format-item )
!
! Not [completely] tested here: format-item-list.
!
! Duplicate statements occur to test unparsing. These duplicates
! should be symantically the same (needs to be checked) and
! should unparse the same way.
!
0010 FORMAT (I12, i1)
0011 FORMAT(I12,I1)
0020 format (1pe12.4)
0021 FORMAT(1P,E12.4)
0030 format (1PE12.4, I10)
0031 FORMAT(1P,E12.4,I10)
0040 format (I12 / 'Dates:', 22I3)
0041 format (I12,/,'Dates:',22I3)
0050 format(//)
0051 format(/,/)
0060 format(:)
0070 format(:::)
0071 format(:,:,:)
0080  format(::)
0081 format(:,:)
0090 format(I12,/,' Dates: ', 2(2I3,I5))
0100 format (I12 / 'Dates:', (2I3, I5))
0101 FORMAT(I12,/,'Dates:', (2I3, I5))
0111 format (I12 / 'Dates:', 2(2I3, I5))
0112 FORMAT(I12,/,'Dates:', 2(2I3,I5))
0121 format (1x, i1, 1x, 'isn''t', 1x, i1)
0122 FORMAT(1X,I1,1X,'isn''t',1X,I1)
0131 format(1x, 2(f10.3, i5))
0132 FORMAT(1X,2(F10.3,I5))
0141 format (1x, f10.3, i5, f10.3, i5, f10.3, i5)
0142 FORMAT(1X,F10.3,I5,F10.3,I5,F10.3,I5)

end

max_stars_repo_name: OpenFortranProject/ofp-sdf

max_issues_repo_name: OpenFortranProject/ofp-sdf

max_forks_repo_name: OpenFortranProject/ofp-sdf


MODULE diahsb
  USE oce
  USE dom_oce
  USE phycst
  USE sbc_oce
  USE sbcrnf
  USE sbcisf
  USE domvvl
  USE traqsr
  USE trabbc
  USE trabbc
  USE restart
  USE bdy_oce, ONLY: ln_bdy
  USE iom
  USE in_out_manager
  USE lib_fortran
  USE lib_mpp
  USE timing
  IMPLICIT NONE
  PRIVATE
  PUBLIC :: dia_hsb
  PUBLIC :: dia_hsb_init
  LOGICAL, PUBLIC :: ln_diahsb
  REAL(KIND = wp) :: surf_tot
  REAL(KIND = wp) :: frc_t, frc_s, frc_v
  REAL(KIND = wp) :: frc_wn_t, frc_wn_s
  REAL(KIND = wp), DIMENSION(:, :), ALLOCATABLE :: surf
  REAL(KIND = wp), DIMENSION(:, :), ALLOCATABLE :: surf_ini, ssh_ini
  REAL(KIND = wp), DIMENSION(:, :), ALLOCATABLE :: ssh_hc_loc_ini, ssh_sc_loc_ini
  REAL(KIND = wp), DIMENSION(:, :, :), ALLOCATABLE :: hc_loc_ini, sc_loc_ini, e3t_ini
  CONTAINS
  SUBROUTINE dia_hsb(kt)
    INTEGER, INTENT(IN) :: kt
    INTEGER :: ji, jj, jk
    REAL(KIND = wp) :: zdiff_hc, zdiff_sc
    REAL(KIND = wp) :: zdiff_hc1, zdiff_sc1
    REAL(KIND = wp) :: zdiff_v1, zdiff_v2
    REAL(KIND = wp) :: zerr_hc1, zerr_sc1
    REAL(KIND = wp) :: zvol_tot
    REAL(KIND = wp) :: z_frc_trd_t, z_frc_trd_s
    REAL(KIND = wp) :: z_frc_trd_v
    REAL(KIND = wp) :: z_wn_trd_t, z_wn_trd_s
    REAL(KIND = wp) :: z_ssh_hc, z_ssh_sc
    REAL(KIND = wp), DIMENSION(jpi, jpj) :: z2d0, z2d1
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpkm1) :: zwrk
    IF (ln_timing) CALL timing_start('dia_hsb')
    !$ACC KERNELS
    tsn(:, :, :, 1) = tsn(:, :, :, 1) * tmask(:, :, :)
    tsb(:, :, :, 1) = tsb(:, :, :, 1) * tmask(:, :, :)
    tsn(:, :, :, 2) = tsn(:, :, :, 2) * tmask(:, :, :)
    tsb(:, :, :, 2) = tsb(:, :, :, 2) * tmask(:, :, :)
    !$ACC END KERNELS
    z_frc_trd_v = r1_rau0 * glob_sum(- (emp(:, :) - rnf(:, :) + fwfisf(:, :)) * surf(:, :))
    z_frc_trd_t = glob_sum(sbc_tsc(:, :, jp_tem) * surf(:, :))
    z_frc_trd_s = glob_sum(sbc_tsc(:, :, jp_sal) * surf(:, :))
    IF (ln_rnf) z_frc_trd_t = z_frc_trd_t + glob_sum(rnf_tsc(:, :, jp_tem) * surf(:, :))
    IF (ln_rnf_sal) z_frc_trd_s = z_frc_trd_s + glob_sum(rnf_tsc(:, :, jp_sal) * surf(:, :))
    IF (ln_isf) z_frc_trd_t = z_frc_trd_t + glob_sum(risf_tsc(:, :, jp_tem) * surf(:, :))
    IF (ln_traqsr) z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * glob_sum(qsr(:, :) * surf(:, :))
    IF (ln_trabbc) z_frc_trd_t = z_frc_trd_t + glob_sum(qgh_trd0(:, :) * surf(:, :))
    IF (ln_linssh) THEN
      IF (ln_isfcav) THEN
        !$ACC KERNELS
        DO ji = 1, jpi
          DO jj = 1, jpj
            z2d0(ji, jj) = surf(ji, jj) * wn(ji, jj, mikt(ji, jj)) * tsb(ji, jj, mikt(ji, jj), jp_tem)
            z2d1(ji, jj) = surf(ji, jj) * wn(ji, jj, mikt(ji, jj)) * tsb(ji, jj, mikt(ji, jj), jp_sal)
          END DO
        END DO
        !$ACC END KERNELS
      ELSE
        !$ACC KERNELS
        z2d0(:, :) = surf(:, :) * wn(:, :, 1) * tsb(:, :, 1, jp_tem)
        z2d1(:, :) = surf(:, :) * wn(:, :, 1) * tsb(:, :, 1, jp_sal)
        !$ACC END KERNELS
      END IF
      z_wn_trd_t = - glob_sum(z2d0)
      z_wn_trd_s = - glob_sum(z2d1)
    END IF
    frc_v = frc_v + z_frc_trd_v * rdt
    frc_t = frc_t + z_frc_trd_t * rdt
    frc_s = frc_s + z_frc_trd_s * rdt
    IF (ln_linssh) THEN
      frc_wn_t = frc_wn_t + z_wn_trd_t * rdt
      frc_wn_s = frc_wn_s + z_wn_trd_s * rdt
    END IF
    zdiff_v1 = glob_sum_full(surf(:, :) * sshn(:, :) - surf_ini(:, :) * ssh_ini(:, :))
    IF (ln_linssh) THEN
      IF (ln_isfcav) THEN
        !$ACC KERNELS
        DO ji = 1, jpi
          DO jj = 1, jpj
            z2d0(ji, jj) = surf(ji, jj) * (tsn(ji, jj, mikt(ji, jj), jp_tem) * sshn(ji, jj) - ssh_hc_loc_ini(ji, jj))
            z2d1(ji, jj) = surf(ji, jj) * (tsn(ji, jj, mikt(ji, jj), jp_sal) * sshn(ji, jj) - ssh_sc_loc_ini(ji, jj))
          END DO
        END DO
        !$ACC END KERNELS
      ELSE
        !$ACC KERNELS
        z2d0(:, :) = surf(:, :) * (tsn(:, :, 1, jp_tem) * sshn(:, :) - ssh_hc_loc_ini(:, :))
        z2d1(:, :) = surf(:, :) * (tsn(:, :, 1, jp_sal) * sshn(:, :) - ssh_sc_loc_ini(:, :))
        !$ACC END KERNELS
      END IF
      z_ssh_hc = glob_sum_full(z2d0)
      z_ssh_sc = glob_sum_full(z2d1)
    END IF
    !$ACC KERNELS
    DO jk = 1, jpkm1
      zwrk(:, :, jk) = (surf(:, :) * e3t_n(:, :, jk) - surf_ini(:, :) * e3t_ini(:, :, jk)) * tmask(:, :, jk)
    END DO
    !$ACC END KERNELS
    zdiff_v2 = glob_sum_full(zwrk(:, :, :))
    !$ACC KERNELS
    DO jk = 1, jpkm1
      zwrk(:, :, jk) = (surf(:, :) * e3t_n(:, :, jk) * tsn(:, :, jk, jp_tem) - surf_ini(:, :) * hc_loc_ini(:, :, jk)) * tmask(:, :, jk)
    END DO
    !$ACC END KERNELS
    zdiff_hc = glob_sum_full(zwrk(:, :, :))
    !$ACC KERNELS
    DO jk = 1, jpkm1
      zwrk(:, :, jk) = (surf(:, :) * e3t_n(:, :, jk) * tsn(:, :, jk, jp_sal) - surf_ini(:, :) * sc_loc_ini(:, :, jk)) * tmask(:, :, jk)
    END DO
    !$ACC END KERNELS
    zdiff_sc = glob_sum_full(zwrk(:, :, :))
    zdiff_v1 = zdiff_v1 - frc_v
    IF (.NOT. ln_linssh) zdiff_v2 = zdiff_v2 - frc_v
    zdiff_hc = zdiff_hc - frc_t
    zdiff_sc = zdiff_sc - frc_s
    IF (ln_linssh) THEN
      zdiff_hc1 = zdiff_hc + z_ssh_hc
      zdiff_sc1 = zdiff_sc + z_ssh_sc
      zerr_hc1 = z_ssh_hc - frc_wn_t
      zerr_sc1 = z_ssh_sc - frc_wn_s
    END IF
    !$ACC KERNELS
    DO jk = 1, jpkm1
      zwrk(:, :, jk) = surf(:, :) * e3t_n(:, :, jk) * tmask(:, :, jk)
    END DO
    !$ACC END KERNELS
    zvol_tot = glob_sum_full(zwrk(:, :, :))
    CALL iom_put('bgfrcvol', frc_v * 1.E-9)
    CALL iom_put('bgfrctem', frc_t * rau0 * rcp * 1.E-20)
    CALL iom_put('bgfrchfx', frc_t * rau0 * rcp / (surf_tot * kt * rdt))
    CALL iom_put('bgfrcsal', frc_s * 1.E-9)
    IF (.NOT. ln_linssh) THEN
      CALL iom_put('bgtemper', zdiff_hc / zvol_tot)
      CALL iom_put('bgsaline', zdiff_sc / zvol_tot)
      CALL iom_put('bgheatco', zdiff_hc * 1.E-20 * rau0 * rcp)
      CALL iom_put('bgheatfx', zdiff_hc * rau0 * rcp / (surf_tot * kt * rdt))
      CALL iom_put('bgsaltco', zdiff_sc * 1.E-9)
      CALL iom_put('bgvolssh', zdiff_v1 * 1.E-9)
      CALL iom_put('bgvole3t', zdiff_v2 * 1.E-9)
      IF (kt == nitend .AND. lwp) THEN
        WRITE(numout, FMT = *)
        WRITE(numout, FMT = *) 'dia_hsb : last time step hsb diagnostics: at it= ', kt, ' date= ', ndastp
        WRITE(numout, FMT = *) '~~~~~~~'
        WRITE(numout, FMT = *) '   Temperature drift = ', zdiff_hc / zvol_tot, ' C'
        WRITE(numout, FMT = *) '   Salinity    drift = ', zdiff_sc / zvol_tot, ' PSU'
        WRITE(numout, FMT = *) '   volume ssh  drift = ', zdiff_v1 * 1.E-9, ' km^3'
        WRITE(numout, FMT = *) '   volume e3t  drift = ', zdiff_v2 * 1.E-9, ' km^3'
      END IF
    ELSE
      CALL iom_put('bgtemper', zdiff_hc1 / zvol_tot)
      CALL iom_put('bgsaline', zdiff_sc1 / zvol_tot)
      CALL iom_put('bgheatco', zdiff_hc1 * 1.E-20 * rau0 * rcp)
      CALL iom_put('bgheatfx', zdiff_hc1 * rau0 * rcp / (surf_tot * kt * rdt))
      CALL iom_put('bgsaltco', zdiff_sc1 * 1.E-9)
      CALL iom_put('bgvolssh', zdiff_v1 * 1.E-9)
      CALL iom_put('bgmistem', zerr_hc1 / zvol_tot)
      CALL iom_put('bgmissal', zerr_sc1 / zvol_tot)
    END IF
    IF (lrst_oce) CALL dia_hsb_rst(kt, 'WRITE')
    IF (ln_timing) CALL timing_stop('dia_hsb')
  END SUBROUTINE dia_hsb
  SUBROUTINE dia_hsb_rst(kt, cdrw)
    INTEGER, INTENT(IN) :: kt
    CHARACTER(LEN = *), INTENT(IN) :: cdrw
    INTEGER :: ji, jj, jk
    IF (TRIM(cdrw) == 'READ') THEN
      IF (ln_rstart) THEN
        IF (lwp) WRITE(numout, FMT = *)
        IF (lwp) WRITE(numout, FMT = *) '   dia_hsb_rst : read hsb restart at it= ', kt, ' date= ', ndastp
        IF (lwp) WRITE(numout, FMT = *)
        CALL iom_get(numror, 'frc_v', frc_v, ldxios = lrxios)
        CALL iom_get(numror, 'frc_t', frc_t, ldxios = lrxios)
        CALL iom_get(numror, 'frc_s', frc_s, ldxios = lrxios)
        IF (ln_linssh) THEN
          CALL iom_get(numror, 'frc_wn_t', frc_wn_t, ldxios = lrxios)
          CALL iom_get(numror, 'frc_wn_s', frc_wn_s, ldxios = lrxios)
        END IF
        CALL iom_get(numror, jpdom_autoglo, 'surf_ini', surf_ini, ldxios = lrxios)
        CALL iom_get(numror, jpdom_autoglo, 'ssh_ini', ssh_ini, ldxios = lrxios)
        CALL iom_get(numror, jpdom_autoglo, 'e3t_ini', e3t_ini, ldxios = lrxios)
        CALL iom_get(numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini, ldxios = lrxios)
        CALL iom_get(numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini, ldxios = lrxios)
        IF (ln_linssh) THEN
          CALL iom_get(numror, jpdom_autoglo, 'ssh_hc_loc_ini', ssh_hc_loc_ini, ldxios = lrxios)
          CALL iom_get(numror, jpdom_autoglo, 'ssh_sc_loc_ini', ssh_sc_loc_ini, ldxios = lrxios)
        END IF
      ELSE
        IF (lwp) WRITE(numout, FMT = *)
        IF (lwp) WRITE(numout, FMT = *) '   dia_hsb_rst : initialise hsb at initial state '
        IF (lwp) WRITE(numout, FMT = *)
        !$ACC KERNELS
        surf_ini(:, :) = e1e2t(:, :) * tmask_i(:, :)
        ssh_ini(:, :) = sshn(:, :)
        DO jk = 1, jpk
          e3t_ini(:, :, jk) = e3t_n(:, :, jk) * tmask(:, :, jk)
          hc_loc_ini(:, :, jk) = tsn(:, :, jk, jp_tem) * e3t_n(:, :, jk) * tmask(:, :, jk)
          sc_loc_ini(:, :, jk) = tsn(:, :, jk, jp_sal) * e3t_n(:, :, jk) * tmask(:, :, jk)
        END DO
        frc_v = 0._wp
        frc_t = 0._wp
        frc_s = 0._wp
        !$ACC END KERNELS
        IF (ln_linssh) THEN
          IF (ln_isfcav) THEN
            !$ACC KERNELS
            DO ji = 1, jpi
              DO jj = 1, jpj
                ssh_hc_loc_ini(ji, jj) = tsn(ji, jj, mikt(ji, jj), jp_tem) * sshn(ji, jj)
                ssh_sc_loc_ini(ji, jj) = tsn(ji, jj, mikt(ji, jj), jp_sal) * sshn(ji, jj)
              END DO
            END DO
            !$ACC END KERNELS
          ELSE
            !$ACC KERNELS
            ssh_hc_loc_ini(:, :) = tsn(:, :, 1, jp_tem) * sshn(:, :)
            ssh_sc_loc_ini(:, :) = tsn(:, :, 1, jp_sal) * sshn(:, :)
            !$ACC END KERNELS
          END IF
          frc_wn_t = 0._wp
          frc_wn_s = 0._wp
        END IF
      END IF
    ELSE IF (TRIM(cdrw) == 'WRITE') THEN
      IF (lwp) WRITE(numout, FMT = *)
      IF (lwp) WRITE(numout, FMT = *) '   dia_hsb_rst : write restart at it= ', kt, ' date= ', ndastp
      IF (lwp) WRITE(numout, FMT = *)
      IF (lwxios) CALL iom_swap(cwxios_context)
      CALL iom_rstput(kt, nitrst, numrow, 'frc_v', frc_v, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'frc_t', frc_t, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'frc_s', frc_s, ldxios = lwxios)
      IF (ln_linssh) THEN
        CALL iom_rstput(kt, nitrst, numrow, 'frc_wn_t', frc_wn_t, ldxios = lwxios)
        CALL iom_rstput(kt, nitrst, numrow, 'frc_wn_s', frc_wn_s, ldxios = lwxios)
      END IF
      CALL iom_rstput(kt, nitrst, numrow, 'surf_ini', surf_ini, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'ssh_ini', ssh_ini, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'e3t_ini', e3t_ini, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini, ldxios = lwxios)
      CALL iom_rstput(kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini, ldxios = lwxios)
      IF (ln_linssh) THEN
        CALL iom_rstput(kt, nitrst, numrow, 'ssh_hc_loc_ini', ssh_hc_loc_ini, ldxios = lwxios)
        CALL iom_rstput(kt, nitrst, numrow, 'ssh_sc_loc_ini', ssh_sc_loc_ini, ldxios = lwxios)
      END IF
      IF (lwxios) CALL iom_swap(cxios_context)
    END IF
  END SUBROUTINE dia_hsb_rst
  SUBROUTINE dia_hsb_init
    INTEGER :: ierror, ios
    NAMELIST /namhsb/ ln_diahsb
    IF (lwp) THEN
      WRITE(numout, FMT = *)
      WRITE(numout, FMT = *) 'dia_hsb_init : heat and salt budgets diagnostics'
      WRITE(numout, FMT = *) '~~~~~~~~~~~~ '
    END IF
    REWIND(UNIT = numnam_ref)
    READ(numnam_ref, namhsb, IOSTAT = ios, ERR = 901)
901 IF (ios /= 0) CALL ctl_nam(ios, 'namhsb in reference namelist', lwp)
    REWIND(UNIT = numnam_cfg)
    READ(numnam_cfg, namhsb, IOSTAT = ios, ERR = 902)
902 IF (ios > 0) CALL ctl_nam(ios, 'namhsb in configuration namelist', lwp)
    IF (lwm) WRITE(numond, namhsb)
    IF (lwp) THEN
      WRITE(numout, FMT = *) '   Namelist  namhsb :'
      WRITE(numout, FMT = *) '      check the heat and salt budgets (T) or not (F)       ln_diahsb = ', ln_diahsb
    END IF
    IF (.NOT. ln_diahsb) RETURN
    IF (lwxios) THEN
      CALL iom_set_rstw_var_active('frc_v')
      CALL iom_set_rstw_var_active('frc_t')
      CALL iom_set_rstw_var_active('frc_s')
      CALL iom_set_rstw_var_active('surf_ini')
      CALL iom_set_rstw_var_active('ssh_ini')
      CALL iom_set_rstw_var_active('e3t_ini')
      CALL iom_set_rstw_var_active('hc_loc_ini')
      CALL iom_set_rstw_var_active('sc_loc_ini')
      IF (ln_linssh) THEN
        CALL iom_set_rstw_var_active('ssh_hc_loc_ini')
        CALL iom_set_rstw_var_active('ssh_sc_loc_ini')
        CALL iom_set_rstw_var_active('frc_wn_t')
        CALL iom_set_rstw_var_active('frc_wn_s')
      END IF
    END IF
    ALLOCATE(hc_loc_ini(jpi, jpj, jpk), sc_loc_ini(jpi, jpj, jpk), surf_ini(jpi, jpj), e3t_ini(jpi, jpj, jpk), surf(jpi, jpj), ssh_ini(jpi, jpj), STAT = ierror)
    IF (ierror > 0) THEN
      CALL ctl_stop('dia_hsb_init: unable to allocate hc_loc_ini')
      RETURN
    END IF
    IF (ln_linssh) ALLOCATE(ssh_hc_loc_ini(jpi, jpj), ssh_sc_loc_ini(jpi, jpj), STAT = ierror)
    IF (ierror > 0) THEN
      CALL ctl_stop('dia_hsb: unable to allocate ssh_hc_loc_ini')
      RETURN
    END IF
    !$ACC KERNELS
    surf(:, :) = e1e2t(:, :) * tmask_i(:, :)
    !$ACC END KERNELS
    surf_tot = glob_sum(surf(:, :))
    IF (ln_bdy) CALL ctl_warn('dia_hsb_init: heat/salt budget does not consider open boundary fluxes')
    CALL dia_hsb_rst(nit000, 'READ')
  END SUBROUTINE dia_hsb_init
END MODULE diahsb


text:
        SUBROUTINE SOFIO ( ISOP, IBLKNM, BUF )
      INCLUDE 'GINOX.COM'
      INCLUDE 'DSIOF.COM'
        COMMON / SOFCOM / NFILES, FILNAM( 10 ), FILSIZ( 10 )
        COMMON / SYS    / BLKSIZ, DIRSIZ, SUPSIZ, AVBLKS, HIBLK
        COMMON / SYSTEM / ISYSBF, IWR
        common / sofdsn / sofdsn(10)
        CHARACTER*4       FILNAM
        CHARACTER*80      DSNAME
        character*80      sofdsn
        INTEGER           FILSIZ, HIBLK, BUF(10)
        IF ( LENSOF( 1 ) .NE. 0 ) GO TO 20
        NUMBLK = 1
        IF ( LENWPB .NE. 0 ) NUMBLK = ISYSBF / LENWPB
        DO 10 K = 1, NFILES
        LENSOF( K ) = 0
        DSNAME      = sofdsn(K)
        CALL DSINQR ( DSNAME, ISTAT, ISIZE)
        IF (ISTAT .EQ. 0) GO TO 10
        LENSOF( K ) = FILSIZ( K )
10      CONTINUE
        LASFIL = 0
20      CONTINUE
        IF ( ISOP .EQ. 7 ) GO TO 200
        IBLK   = IBLKNM
        IF ( IBLK .LE. 0 ) GO TO 700
        IFILE  = 0
        DO 50 K = 1, NFILES
        IF ( IBLK .LE. FILSIZ( K ) )GO TO 30
        IBLK   = IBLK - FILSIZ( K )
        GO TO 50
30      IFILE  = K
        GO TO 100
50      CONTINUE
        WRITE( IWR, 9910 ) IBLKNM
9910    FORMAT(' *** SUBSTRUCTURING ERROR - BLOCK NUMBER OUT OF RANGE ',
     *         ' -  BLOCK NUMBER IS ',I10)
        WRITE( IWR, 9915 )
9915    FORMAT( //,' THE FOLLOWING SOF FILES WERE AVAILABLE',//)
        DO 60 K = 1, NFILES
        WRITE( IWR, 9920 ) sofdsn( K ), FILSIZ( K ), LENSOF( K )
9920    FORMAT(' FILE ',A72' HAS ',I10, ' BLOCKS - BLOCKS USED ',I10)
60      CONTINUE
        CALL MESAGE (-61, 0, 0)
100     IF ( LASFIL .EQ. IFILE ) GO TO 120
        IF ( LASFIL .NE. 0 ) CALL DSCLOS ( 90 )
        IALLOC = NUMBLK * FILSIZ(IFILE)
        dsname = sofdsn( IFILE )
        IOP = 0
        CALL DSOPEN (  DSNAME, 90, IOP )
        LASFIL = IFILE
120     IF ( ISOP .EQ. 1 ) GO TO 140
        IF ( ( IBLK - LENSOF( IFILE ) ) .LE. 1 ) GO TO 130
        NUM  = IBLK - LENSOF( IFILE ) - 1
        IF ( NUM .EQ. 0 ) GO TO 130
        DO 125 K = 1, NUM
        LENSOF( IFILE ) = LENSOF( IFILE ) + 1
        CALL DSWRIT ( 90, BUF(4), NBUFF, LENSOF( IFILE ),ICERR )
        IF ( ICERR .NE. 0 ) GO TO 701
125     CONTINUE
130     CONTINUE
        CALL DSWRIT ( 90, BUF(4), NBUFF, IBLK, ICERR )
        IF ( ICERR .NE. 0 ) GO TO 701
        IF ( IBLK   .GT. LENSOF( IFILE ) ) LENSOF( IFILE ) = IBLK
        IF ( IBLKNM .GT. HIBLK ) HIBLK = IBLKNM
        GO TO 700
140     CALL DSREAD ( 90, BUF(4), NBUFF, IBLK )
        GO TO 700
200     CALL DSCLOS( 90 )
        LASFIL = 0
700     GO TO 7000
701     IF ( ICERR .EQ. 28 ) WRITE ( IWR, 901 )
        IF ( ICERR .NE. 28 ) WRITE ( IWR, 902 )
        CALL MESAGE (-61, 0, 0)
901     FORMAT(///,' INSUFFICIENT SPACE FOR SOF FILE ON DEFAULT',
     &             ' DEVICE---JOB ABORTED.')
902     FORMAT(///,' I/O ERROR OCCURRED ON SOF FILE, JOB ABORTED')
7000    RETURN
        END


text:
!------------------------------------------------------------------------------------------
! File:   HMC_Module_Phys_ET.f90
! Author: Fabio Delogu
!
! Created on April 2, 2014, 5:19 PM
!------------------------------------------------------------------------------------------

!------------------------------------------------------------------------------------------
! Module Header
module HMC_Module_Phys_ET

    !------------------------------------------------------------------------------------------
    ! External module(s) for all subroutine in this module
    use HMC_Module_Vars_Loader, only: oHMC_Vars

    use HMC_Module_Tools_Debug

    ! Implicit none for all subroutines in this module
    implicit none
    !------------------------------------------------------------------------------------------

contains

    !------------------------------------------------------------------------------------------
    ! Subroutine to calculate evapotranspiration
    subroutine HMC_Phys_ET_Cpl(iID, iRows, iCols, iTime, iDaySteps, sTime, iNLake, iNDam)

	!------------------------------------------------------------------------------------------
        ! Variable(s) declaration
        integer(kind = 4)           :: iID
        integer(kind = 4)           :: iRows, iCols
        integer(kind = 4)           :: iNLake, iNDam
        integer(kind = 4)           :: iTime, iDaySteps, iStep
        integer(kind = 4)           :: iI, iJ, iL, iD

        real(kind = 4), dimension (iRows, iCols)         :: a2dVarET, a2dVarVTot, a2dVarETPot
        real(kind = 4), dimension (iRows, iCols)         :: a2dVarVRet, a2dVarVTotWP

        real(kind = 4)              :: dVarET, dVarAE, dVarETLake, dVarETTot, dVarETPot

        integer(kind = 4), dimension (iRows, iCols)         :: a2iVarMask, a2iVarChoice
        real(kind = 4), dimension (iRows, iCols)            :: a2dVarDEM, a2dVarAreaCell, a2dVarCtWP

        real(kind = 4), dimension (iRows, iCols)            :: a2dVarAE, a2dVarAEres

        integer(kind = 4), dimension (iNDam, 2)             :: a2iVarXYDam
        real(kind = 4), dimension (iNDam)                   :: a1dVarCodeDam
        real(kind = 4), dimension (iNDam)                   :: a1dVarVDam
        integer(kind = 4), dimension (iNLake, 2)            :: a2iVarXYLake
        real(kind = 4), dimension (iNLake)                  :: a1dVarCodeLake
        real(kind = 4), dimension (iNLake)                  :: a1dVarVLake

        character(len = 19)                                 :: sTime

        character(len = 10)                                 :: sVarAE, sVarET, sVarETTot, sVarETPot
        character(len = 10), parameter                      :: sFMTVarET = "(F7.4)"
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Initialization variable(s)
        a2iVarMask = 0; a2iVarChoice = 0; a2dVarDEM = 0.0; a2dVarAreaCell = 0.0; a2dVarCtWP = 0.0;
        a2dVarAE = 0.0; a2dVarET = 0.0; a2dVarETPot = 0.0;
        a2dVarVRet = 0.0; a2dVarVTot = 0.0; a2dVarVTotWP = 0.0;
        dVarAE = 0.0; dVarET = 0.0; dVarETLake = 0.0; a2dVarAEres = 0.0;
        dVarETTot = 0.0; dVarETPot = 0.0;

        a2iVarXYDam = 0; a1dVarCodeDam = 0.0; a1dVarVDam = 0.0;
        a2iVarXYLake = 0; a1dVarCodeLake = 0.0; a1dVarVLake = 0.0;
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Data static definition
        a2dVarDEM = oHMC_Vars(iID)%a2dDem
        a2iVarMask = oHMC_Vars(iID)%a2iMask
        a2iVarChoice = oHMC_Vars(iID)%a2iChoice
        a2dVarAreaCell = oHMC_Vars(iID)%a2dAreaCell
        a2dVarCtWP = oHMC_Vars(iID)%a2dCtWP

        ! Lake(s) and dam(s) variable(s)
        a2iVarXYDam = oHMC_Vars(iID)%a2iXYDam
        a1dVarCodeDam = oHMC_Vars(iID)%a1dCodeDam
        a1dVarVDam = oHMC_Vars(iID)%a1dVDam
        a2iVarXYLake = oHMC_Vars(iID)%a2iXYLake
        a1dVarCodeLake = oHMC_Vars(iID)%a1dCodeLake
        a1dVarVLake = oHMC_Vars(iID)%a1dVLake

        ! Extracting dynamic variable(s)
        a2dVarVRet = oHMC_Vars(iID)%a2dVRet

        ! Extracting dynamic state variable(s)
        a2dVarVTot = oHMC_Vars(iID)%a2dVTot         ! Total soil volume
        a2dVarET = oHMC_Vars(iID)%a2dET             ! Evapotranspiration ==> from LST phys or forcing dataset
        a2dVarETPot = oHMC_Vars(iID)%a2dETPot       ! Potential Evapotranspiration == from LST phys or forcing dataset

        ! Compute soil water content at wilting point - minimum threshold for evapotranspiration
        a2dVarVTotWP = oHMC_Vars(iID)%a2dVTotWP

        ! Extracting checking variable(s)
        dVarETTot = oHMC_Vars(iID)%dVarETTot

        ! Info start
        call mprintf(.true., iINFO_Verbose, ' Phys :: Evapotranspiration ... ' )
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Debug
        if (iDEBUG.gt.0) then
            call mprintf(.true., iINFO_Extra, ' ========= EVAPOTRANSPIRATION START =========== ')
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarET, a2iVarMask, 'ET START ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarAE, a2iVarMask, 'AE START ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarETPot, a2iVarMask, 'ET.POT START ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarVTot, a2iVarMask, 'VTOT START ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRet, a2iVarMask, 'VRET START ') )
            call mprintf(.true., iINFO_Extra, '')
        endif
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Calculating control variable(s)
        dVarAE = sum(a2dVarET, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarET = sum(a2dVarET, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarETPot = sum(a2dVarETPot, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarETTot = dVarETTot + dVarAE

        ! ET information time step
        write(sVarAE, sFMTVarET) dVarAE
        write(sVarETPot, sFMTVarET) dVarETPot
        write(sVarET, sFMTVarET) dVarET
        write(sVarETTot, sFMTVarET) dVarETTot
        call mprintf(.true., iINFO_Basic, ' Phys :: EVT START :: AvgValue :: '// &
                                          ' AEvt: '//sVarAE//' [mm] '// &
                                          ' PEvt: '//sVarETPot//' [mm] '// &
                                          ' Evt: '//sVarET//' [mm] '// &
                                          ' Evt Tot: '//sVarETTot//' [mm]')
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Lake and dam ET updating

        ! Dam updating
        if (iNDam .gt. 0) then

            ! Cycle on dam(s)
            do iD = 1, iNDam
                iI = 0; iJ = 0;
                iI = a2iVarXYDam(iD, 2)
		iJ = a2iVarXYDam(iD, 1)

                ! Compute dam volume
                if (a1dVarCodeDam(iD) .gt. 0) then
                    ! Distributed dam-lake
                    dVarETLake = 0.0
                    dVarETLake = sum(sum( a2dVarETPot, DIM = 1, mask=a2iVarChoice .eq. a1dVarCodeDam(iD)))
                    ! Update dam volume
                    a1dVarVDam(iD) = a1dVarVDam(iD) - dVarETLake/1000*a2dVarAreaCell(iI, iJ) !in m^3
                endif

                if (a1dVarVDam(iD) .lt. 0.0) a1dVarVDam(iD) = 0.0
                dVarETLake = 0.0

            enddo
        endif

        ! Lake updating
        if (iNLake .gt. 0 ) then

            ! Cycle on lake(s)
            do iL = 1, iNLake

                iI = 0; iJ = 0;
                iI = a2iVarXYLake(iL,2)
                iJ = a2iVarXYLake(iL,1)

                ! Compute lake volume
                if (a1dVarCodeLake(iL) .gt. 0) then
                    ! Distributed lake
                    dVarETLake = 0.0
                    dVarETLake = sum(sum( a2dVarETPot, DIM = 1, mask=a2iVarChoice .eq. a1dVarCodeLake(iL)))
                    ! Update dam volume
                    a1dVarVLake(iL) = a1dVarVLake(iL) - dVarETLake/1000*a2dVarAreaCell(iI, iJ) !in m^3
                endif

                if (a1dVarVLake(iL) .lt. 0.0) a1dVarVLake(iL) = 0.0
                dVarETLake = 0.0

            enddo
        endif
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Passing evapotranspiration (ET) to Actual evapotranspiration (AE)
        where (a2dVarDEM.gt.0.0)
            a2dVarAE = a2dVarET
        endwhere
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Calculating retention volume
        where( (a2dVarVRet.gt.0.0) .and. (a2dVarVRet.ge.a2dVarETPot) .and.  &
               (a2dVarDEM.gt.0.0) .and. (a2iVarChoice.le.1) )

            a2dVarVRet = a2dVarVRet - a2dVarETPot
            a2dVarAE = a2dVarETPot

        elsewhere( (a2dVarVRet.gt.0.0) .and. (a2dVarVRet.lt.a2dVarETPot) .and. &
                (a2dVarDEM.gt.0.0) .and. (a2iVarChoice.le.1) )

            ! Compute residual evapotranspiration demand
            a2dVarAEres = a2dVarETPot - a2dVarVRet

            where( a2dVarAEres .lt. a2dVarAE)
                a2dVarAE = a2dVarAEres
            endwhere

            where( a2dVarAE .lt. 0.0)
                a2dVarAE = 0.0
            endwhere

            where( (a2dVarVTot - a2dVarVTotWP) .gt. a2dVarAE )
                a2dVarVTot = a2dVarVTot - a2dVarAE
                a2dVarAE = a2dVarVRet + a2dVarAE
            elsewhere ( (a2dVarVTot - a2dVarVTotWP) .le. a2dVarAE .and. (a2dVarVTot - a2dVarVTotWP) .gt. 0.0 )
                a2dVarAE = a2dVarVRet + (a2dVarVTot - a2dVarVTotWP)
                a2dVarVTot = a2dVarVTotWP

            elsewhere ( (a2dVarVTot-a2dVarVTotWP) .le. 0.0 ) ! to account also for VTot<Vwp situations
                a2dVarAE = 0.0
            endwhere
            a2dVarVRet = 0.0

        elsewhere (a2dVarDEM.gt.0.0 .and. a2iVarChoice.le.1) ! Retention == 0.0 not on lakes

            where( a2dVarAE .lt. 0.0)
                a2dVarAE = 0.0
            endwhere

            where ((a2dVarVTot - a2dVarVTotWP).gt.a2dVarAE)
                ! tolgo evt da a2dV solo quando "non piove" cio� a2dRetention=0
                ! quando piove l'evaporazione dal suolo � trascurabile
                a2dVarVTot = a2dVarVTot - a2dVarAE
            elsewhere ((a2dVarVTot-a2dVarVTotWP).le.a2dVarAE .and. (a2dVarVTot-a2dVarVTotWP) .gt. 0.0 )
                a2dVarAE = a2dVarVTot - a2dVarVTotWP
                a2dVarVTot = a2dVarVTotWP
            elsewhere ( (a2dVarVTot-a2dVarVTotWP) .le. 0.0 ) ! to account also for VTot<Vwp situations
                a2dVarAE = 0.0
            endwhere

        endwhere
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Calculating control variable(s)
        dVarAE = sum(a2dVarAE, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarET = sum(a2dVarET, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarETPot = sum(a2dVarETPot, mask=a2dVarDem.gt.0.0)/max(1,count(a2dVarDem.gt.0.0))
        dVarETTot = dVarETTot + dVarAE

        ! ET information time step
        write(sVarAE, sFMTVarET) dVarAE
        write(sVarETPot, sFMTVarET) dVarETPot
        write(sVarET, sFMTVarET) dVarET
        write(sVarETTot, sFMTVarET) dVarETTot
        call mprintf(.true., iINFO_Basic, ' Phys :: EVT END :: AvgValue :: '// &
                                          ' AEvt: '//sVarAE//' [mm] '// &
                                          ' PEvt: '//sVarETPot//' [mm] '// &
                                          ' Evt: '//sVarET//' [mm] '// &
                                          ' Evt Tot: '//sVarETTot//' [mm]')
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Debug
        if (iDEBUG.gt.0) then
            call mprintf(.true., iINFO_Extra, '')
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarET, a2iVarMask, 'ET END ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarAE, a2iVarMask, 'AE END ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarETPot, a2iVarMask, 'ET.POT END ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarVTot, a2iVarMask, 'VTOT END ') )
            call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRet, a2iVarMask, 'VRET END ') )
            call mprintf(.true., iINFO_Extra, ' ========= EVAPOTRANSPIRATION END =========== ')
        endif
        !------------------------------------------------------------------------------------------

        !------------------------------------------------------------------------------------------
        ! Copy instantaneous AE & AEpot to 3D var
        ! AE - 3D
        if (all(oHMC_Vars(iID)%a3dAE.lt.0.0))then
            oHMC_Vars(iID)%a3dAE(:,:,int(iDaySteps)) =  a2dVarAE
        else
            ! Re-initializing
            do iStep=2, int(iDaySteps)
                oHMC_Vars(iID)%a3dAE(:,:,int(iStep-1)) = oHMC_Vars(iID)%a3dAE(:,:,int(iStep))
            enddo
            ! Updating with new field
            where(oHMC_Vars(iID)%a2dDEM.gt.0.0)
                oHMC_Vars(iID)%a3dAE(:,:,int(iDaySteps)) =  a2dVarAE
            elsewhere
                oHMC_Vars(iID)%a3dAE(:,:,int(iDaySteps)) = -9999.0
            endwhere
        endif

        ! AEpot - 3D
        if (all(oHMC_Vars(iID)%a3dAEpot.lt.0.0))then
            oHMC_Vars(iID)%a3dAEpot(:,:,int(iDaySteps)) =  a2dVarETPot
        else
            ! Re-initializing
            do iStep=2, int(iDaySteps)
                oHMC_Vars(iID)%a3dAEpot(:,:,int(iStep-1)) = oHMC_Vars(iID)%a3dAEpot(:,:,int(iStep))
            enddo
            ! Updating with new field
            where(oHMC_Vars(iID)%a2dDEM.gt.0.0)
                oHMC_Vars(iID)%a3dAEpot(:,:,int(iDaySteps)) =  a2dVarETPot
            elsewhere
                oHMC_Vars(iID)%a3dAEpot(:,:,int(iDaySteps)) = -9999.0
            endwhere
        endif

        ! Passing variable(s) to global declaration
        oHMC_Vars(iID)%a2dAE = a2dVarAE
        oHMC_Vars(iID)%a2dAEPot = a2dVarETPot
        oHMC_Vars(iID)%a2dVRet = a2dVarVRet
        oHMC_Vars(iID)%a1dVLake = a1dVarVLake
        oHMC_Vars(iID)%a1dVDam = a1dVarVDam

        ! Updating state variable(s): ET and Total Volume
        oHMC_Vars(iID)%a2dET = 0.0              ! Re-initializing ET
        oHMC_Vars(iID)%a2dETPot = 0.0           ! Re-initializing ET pot
        oHMC_Vars(iID)%a2dVTot = a2dVarVTot     ! Updating total volume

        oHMC_Vars(iID)%dVarETTot = dVarETTot    ! Check mean cumulated ET

        ! Info end
        call mprintf(.true., iINFO_Verbose, ' Phys :: Evapotranspiration ... OK ' )
        !------------------------------------------------------------------------------------------

    end subroutine HMC_Phys_ET_Cpl
    !------------------------------------------------------------------------------------------

end module HMC_Module_Phys_ET
!------------------------------------------------------------------------------------------

max_stars_repo_name: c-hydro/hmc-dev

max_issues_repo_name: c-hydro/hmc-dev

max_forks_repo_name: c-hydro/hmc-dev