doomsday.f
1 ! @section Synopsis
2 !
3 ! Doomsday algorithm giving the day of the week for a date.
4 !
5 ! @author J. Marcel van der Veer
6 !
7 ! @section copyright
8 !
9 ! This file is part of VIF - vintage fortran compiler.
10 ! Copyright 2020-2025 J. Marcel van der Veer <algol68g@xs4all.nl>.
11 !
12 ! @section license
13 !
14 ! This program is free software; you can redistribute it and/or modify it
15 ! under the terms of the gnu general public license as published by the
16 ! free software foundation; either version 3 of the license, or
17 ! (at your option) any later version.
18 !
19 ! This program is distributed in the hope that it will be useful, but
20 ! without any warranty; without even the implied warranty of merchantability
21 ! or fitness for a particular purpose. See the GNU general public license for
22 ! more details. You should have received a copy of the GNU general public
23 ! license along with this program. If not, see <http://www.gnu.org/licenses/>.
24 !
25 SUBROUTINE doomsday_gregorian (y, m, d, w)
26
27 c DOOMSDAY_GREGORIAN(): weekday given any date in Gregorian calendar.
28 c
29 c Discussion:
30 c
31 c This procedure does not include any procedure to switch to the Julian
32 c calendar for dates early enough that that calendar was used instead.
33 c
34 c Licensing:
35 c
36 c This code is distributed under the MIT license.
37 c
38 c Modified:
39 c
40 c 28 May 2012
41 c
42 c Author:
43 c
44 c John Burkardt
45 c
46 c Reference:
47 c
48 c John Conway,
49 c Tomorrow is the Day After Doomsday,
50 c Eureka,
51 c Volume 36, October 1973, pages 28-31.
52 c
53 c Parameters:
54 c
55 c Input, integer Y, M, D, the year, month and day of the date.
56 c Note that the year must be positive.
57 c
58 c Output, integer W, the weekday of the date.
59 c
60 IMPLICIT none
61
62 INTEGER anchor(4)
63 INTEGER c
64 INTEGER d
65 INTEGER drd
66 INTEGER drdr
67 INTEGER i4_wrap
68 INTEGER l
69 INTEGER m
70 INTEGER mdoom(12)
71 INTEGER w
72 INTEGER y
73 INTEGER ydoom
74 LOGICAL year_is_leap_gregorian
75 INTEGER yy
76 INTEGER yy12d
77 INTEGER yy12r
78 INTEGER yy12r4d
79
80 SAVE anchor
81 SAVE mdoom
82
83 DATA anchor / 1, 6, 4, 3 /
84 DATA mdoom / 3, 28, 0, 4, 9, 6, 11, 8, 5, 10, 7, 12 /
85 c
86 c Refuse to handle Y <= 0.
87 c
88 IF (y .le. 0) THEN
89 write (*, '(a)') ' '
90 write (*, '(a)') 'DOOMSDAY_GREGORIAN - Fatal error!'
91 write (*, '(a)') ' Y <= 0.'
92 STOP
93 END IF
94 c
95 c Determine the century C.
96 c
97 c = y / 100
98 c
99 c Determine the last two digits of the year, YY
100 c
101 yy = mod (y, 100)
102 c
103 c Divide the last two digits of the year by 12.
104 c
105 yy12d = yy / 12
106 yy12r = mod (yy, 12)
107 yy12r4d = yy12r / 4
108 drd = yy12d + yy12r + yy12r4d
109
110 drdr = mod (drd, 7)
111 ydoom = anchor(mod (c-1, 4) + 1) + drdr
112 ydoom = i4_wrap (ydoom, 1, 7)
113 c
114 c If M = 1 or 2, and Y is a leap year, add 1.
115 c
116 IF ((m .eq. 1 .or. m .eq. 2) .and.
117 & year_is_leap_gregorian (y)) THEN
118 l = 1
119 ELSE
120 l = 0
121 END IF
122
123 w = ydoom + (d - mdoom(m) - l)
124 w = i4_wrap (w, 1, 7)
125
126
127 RETURN
128 END
129
130 FUNCTION i4_modp (i, j)
131
132 cc I4_MODP returns the nonnegative remainder of integer division.
133 c
134 c Discussion:
135 c
136 c If
137 c NREM = I4_MODP (I, J)
138 c NMULT = (I - NREM) / J
139 c then
140 c I = J * NMULT + NREM
141 c where NREM is always nonnegative.
142 c
143 c The MOD function computes a result with the same sign as the
144 c quantity being divided. Thus, suppose you had an angle A,
145 c and you wanted to ensure that it was between 0 and 360.
146 c Then mod(A,360) would do, if A was positive, but if A
147 c was negative, your result would be between -360 and 0.
148 c
149 c On the other hand, I4_MODP(A,360) is between 0 and 360, always.
150 c
151 c Example:
152 c
153 c I J MOD I4_MODP Factorization
154 c
155 c 107 50 7 7 107 = 2 * 50 + 7
156 c 107 -50 7 7 107 = -2 * -50 + 7
157 c -107 50 -7 43 -107 = -3 * 50 + 43
158 c -107 -50 -7 43 -107 = 3 * -50 + 43
159 c
160 c Licensing:
161 c
162 c This code is distributed under the MIT license.
163 c
164 c Modified:
165 c
166 c 30 December 2006
167 c
168 c Author:
169 c
170 c John Burkardt
171 c
172 c Parameters:
173 c
174 c Input, integer I, the number to be divided.
175 c
176 c Input, integer J, the number that divides I.
177 c
178 c Output, integer I4_MODP, the nonnegative remainder when I is
179 c divided by J.
180 c
181 IMPLICIT none
182
183 INTEGER i
184 INTEGER i4_modp
185 INTEGER j
186 INTEGER value
187
188 IF (j .eq. 0) THEN
189 CALL xerabt ('i4_modp: illegal divisor', 1)
190 END IF
191
192 value = mod (i, j)
193
194 IF (value .lt. 0) THEN
195 value = value + abs (j)
196 END IF
197
198 i4_modp = value
199
200 RETURN
201 END
202
203 FUNCTION i4_wrap (ival, ilo, ihi)
204
205 c I4_WRAP forces an I4 to lie between given limits by wrapping.
206 c
207 c Example:
208 c
209 c ILO = 4, IHI = 8
210 c
211 c I Value
212 c
213 c -2 8
214 c -1 4
215 c 0 5
216 c 1 6
217 c 2 7
218 c 3 8
219 c 4 4
220 c 5 5
221 c 6 6
222 c 7 7
223 c 8 8
224 c 9 4
225 c 10 5
226 c 11 6
227 c 12 7
228 c 13 8
229 c 14 4
230 c
231 c Licensing:
232 c
233 c This code is distributed under the MIT license.
234 c
235 c Modified:
236 c
237 c 31 December 2006
238 c
239 c Author:
240 c
241 c John Burkardt
242 c
243 c Parameters:
244 c
245 c Input, integer IVAL, an integer value.
246 c
247 c Input, integer ILO, IHI, the desired bounds for the integer value.
248 c
249 c Output, integer I4_WRAP, a "wrapped" version of IVAL.
250 c
251 IMPLICIT none
252
253 INTEGER i4_modp
254 INTEGER i4_wrap
255 INTEGER ihi
256 INTEGER ilo
257 INTEGER ival
258 INTEGER jhi
259 INTEGER jlo
260 INTEGER value
261 INTEGER wide
262
263 jlo = min (ilo, ihi)
264 jhi = max (ilo, ihi)
265
266 wide = jhi - jlo + 1
267
268 IF (wide .eq. 1) THEN
269 value = jlo
270 ELSE
271 value = jlo + i4_modp (ival - jlo, wide)
272 END IF
273
274 i4_wrap = value
275
276 RETURN
277 END
278
279 SUBROUTINE weekday_to_name_common (w, s)
280
281 c WEEKDAY_TO_NAME_COMMON returns the name of a Common weekday.
282 c
283 c Licensing:
284 c
285 c This code is distributed under the MIT license.
286 c
287 c Modified:
288 c
289 c 28 August 1999
290 c
291 c Author:
292 c
293 c John Burkardt
294 c
295 c Parameters:
296 c
297 c Input, integer W, the weekday index.
298 c
299 c Output, character * (*) S, the weekday name.
300 c
301 IMPLICIT none
302
303 CHARACTER * 9 name(7)
304 CHARACTER * (*) s
305 INTEGER w
306 INTEGER w2
307
308 SAVE name
309
310 DATA name /
311 & 'Sunday ', 'Monday ', 'Tuesday ', 'Wednesday',
312 & 'Thursday ', 'Friday ', 'Saturday ' /
313 c
314 c Make a local copy of the weekday number.
315 c
316 w2 = w
317 c
318 c Return the weekday name.
319 c
320 s = name (w2)
321
322 RETURN
323 END
324
325 SUBROUTINE weekday_values (n_data, y, m, d, w)
326
327 c WEEKDAY_VALUES returns the day of the week for various dates.
328 c
329 c Discussion:
330 c
331 c The CE or Common Era calendar is used, under the
332 c hybrid Julian/Gregorian Calendar, with a transition from Julian
333 c to Gregorian. The day after 04 October 1582 was 15 October 1582.
334 c
335 c The year before 1 AD or CE is 1 BC or BCE. In this data set,
336 c years BC/BCE are indicated by a negative year value.
337 c
338 c Licensing:
339 c
340 c This code is distributed under the MIT license.
341 c
342 c Modified:
343 c
344 c 26 May 2012
345 c
346 c Author:
347 c
348 c John Burkardt
349 c
350 c Reference:
351 c
352 c Edward Reingold, Nachum Dershowitz,
353 c Calendrical Calculations: The Millennium Edition,
354 c Cambridge University Press, 2001,
355 c ISBN: 0 521 77752 6
356 c LC: CE12.R45.
357 c
358 c Parameters:
359 c
360 c Input/output, integer N_DATA. The user sets N_DATA to 0
361 c before the first call. On each call, the routine increments N_DATA by 1,
362 c and returns the corresponding data; when there is no more data, the
363 c output value of N_DATA will be 0 again.
364 c
365 c Output, integer Y, M, D, the Common Era date.
366 c
367 c Output, integer W, the day of the week. Sunday = 1.
368 c
369 IMPLICIT none
370
371 INTEGER n_max
372 PARAMETER (n_max = 34)
373
374 INTEGER d
375 INTEGER d_vec(n_max)
376 INTEGER m
377 INTEGER m_vec(n_max)
378 INTEGER n_data
379 INTEGER w
380 INTEGER w_vec(n_max)
381 INTEGER y
382 INTEGER y_vec(n_max)
383
384 SAVE d_vec
385 SAVE m_vec
386 SAVE w_vec
387 SAVE y_vec
388
389 DATA d_vec /
390 & 30,
391 & 8,
392 & 26,
393 & 3,
394 & 7,
395 & 18,
396 & 7,
397 & 19,
398 & 14,
399 & 18,
400 & 16,
401 & 3,
402 & 26,
403 & 20,
404 & 4,
405 & 25,
406 & 31,
407 & 9,
408 & 24,
409 & 10,
410 & 30,
411 & 24,
412 & 19,
413 & 2,
414 & 27,
415 & 19,
416 & 25,
417 & 29,
418 & 19,
419 & 7,
420 & 17,
421 & 25,
422 & 10,
423 & 18 /
424 DATA m_vec /
425 & 7,
426 & 12,
427 & 9,
428 & 10,
429 & 1,
430 & 5,
431 & 11,
432 & 4,
433 & 10,
434 & 5,
435 & 3,
436 & 3,
437 & 3,
438 & 4,
439 & 6,
440 & 1,
441 & 3,
442 & 9,
443 & 2,
444 & 6,
445 & 6,
446 & 7,
447 & 6,
448 & 8,
449 & 3,
450 & 4,
451 & 8,
452 & 9,
453 & 4,
454 & 10,
455 & 3,
456 & 2,
457 & 11,
458 & 7 /
459 DATA w_vec /
460 & 1,
461 & 4,
462 & 4,
463 & 1,
464 & 4,
465 & 2,
466 & 7,
467 & 1,
468 & 7,
469 & 1,
470 & 6,
471 & 7,
472 & 6,
473 & 1,
474 & 1,
475 & 4,
476 & 7,
477 & 7,
478 & 7,
479 & 4,
480 & 1,
481 & 6,
482 & 1,
483 & 2,
484 & 4,
485 & 1,
486 & 1,
487 & 2,
488 & 2,
489 & 5,
490 & 3,
491 & 1,
492 & 4,
493 & 1 /
494 DATA y_vec /
495 & - 587,
496 & - 169,
497 & 70,
498 & 135,
499 & 470,
500 & 576,
501 & 694,
502 & 1013,
503 & 1066,
504 & 1096,
505 & 1190,
506 & 1240,
507 & 1288,
508 & 1298,
509 & 1391,
510 & 1436,
511 & 1492,
512 & 1553,
513 & 1560,
514 & 1648,
515 & 1680,
516 & 1716,
517 & 1768,
518 & 1819,
519 & 1839,
520 & 1903,
521 & 1929,
522 & 1941,
523 & 1943,
524 & 1943,
525 & 1992,
526 & 1996,
527 & 2038,
528 & 2094 /
529
530 IF (n_data .lt. 0) THEN
531 n_data = 0
532 END IF
533
534 n_data = n_data + 1
535
536 IF (n_max .lt. n_data) THEN
537 n_data = 0
538 y = 0
539 m = 0
540 d = 0
541 w = 0
542 ELSE
543 y = y_vec(n_data)
544 m = m_vec(n_data)
545 d = d_vec(n_data)
546 w = w_vec(n_data)
547 END IF
548
549 RETURN
550 END
551
552 FUNCTION year_is_leap_gregorian (y)
553
554 c YEAR_IS_LEAP_GREGORIAN returns TRUE if the Gregorian year was a leap year.
555 c
556 c Licensing:
557 c
558 c This code is distributed under the MIT license.
559 c
560 c Modified:
561 c
562 c 21 May 2000
563 c
564 c Author:
565 c
566 c John Burkardt
567 c
568 c Parameters:
569 c
570 c Input, integer Y, the year to be checked.
571 c
572 c Output, logical YEAR_IS_LEAP_GREGORIAN, TRUE if the year was a leap year,
573 c FALSE otherwise.
574 c
575 IMPLICIT none
576
577 INTEGER y
578 LOGICAL year_is_leap_gregorian
579
580 IF (y .le. 0) THEN
581 CALL xerabt ('year_is_leap_gregorian: invalid year', 1)
582 END IF
583
584 IF (mod (y, 400) .eq. 0) THEN
585 year_is_leap_gregorian = .true.
586 ELSE IF (mod (y, 100) .eq. 0) then
587 year_is_leap_gregorian = .false.
588 ELSE IF (mod (y, 4) .eq. 0) then
589 year_is_leap_gregorian = .true.
590 ELSE
591 year_is_leap_gregorian = .false.
592 END IF
593
594 RETURN
595 END
|
© 2002-2025 J.M. van der Veer (jmvdveer@xs4all.nl)
|
