1 /* -- translated by f2c (version 20100827).
2 You must link the resulting object file with libf2c:
3 on Microsoft Windows system, link with libf2c.lib;
4 on Linux or Unix systems, link with .../path/to/libf2c.a -lm
5 or, if you install libf2c.a in a standard place, with -lf2c -lm
6 -- in that order, at the end of the command line, as in
7 cc *.o -lf2c -lm
8 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
9
10 http://www.netlib.org/f2c/libf2c.zip
11 */
12
13 #include "libtinyf2c.h"
14
15 /* *********************************************************************** */
16 /* INITIALIZE MATRICES FOR INTERPOLATION */
17 /* *********************************************************************** */
18 /* Copyright (C) 2000-2004 GFD Dennou Club. All rights reserved. */
19 /* ----------------------------------------------------------------------- */
shiniz_(integer * jm,real * x,real * y,real * z__)20 /* Subroutine */ int shiniz_(integer *jm, real *x, real *y, real *z__)
21 {
22 /* System generated locals */
23 integer z_dim1, z_dim2, z_offset, i__1, i__2;
24
25 /* Builtin functions */
26 double sin(doublereal), tan(doublereal);
27
28 /* Local variables */
29 static integer j1, j2;
30 static real t1, t2;
31
32 /* Parameter adjustments */
33 z_dim1 = *jm;
34 z_dim2 = *jm - 0 + 1;
35 z_offset = 1 + z_dim1 * (0 + z_dim2);
36 z__ -= z_offset;
37
38 /* Function Body */
39 i__1 = *jm;
40 for (j1 = 1; j1 <= i__1; ++j1) {
41 i__2 = *jm;
42 for (j2 = 0; j2 <= i__2; ++j2) {
43 t1 = x[j1] - y[j2];
44 t2 = x[j1] + y[j2];
45 z__[j1 + (j2 + z_dim2) * z_dim1] = (sin((*jm << 1) * t1) / tan(t1)
46 + sin((*jm << 1) * t2) / tan(t2)) / (*jm << 1);
47 /* L10: */
48 }
49 z__[j1 + z_dim2 * z_dim1] /= 2;
50 z__[j1 + (*jm + z_dim2) * z_dim1] /= 2;
51 /* L11: */
52 }
53 i__1 = *jm;
54 for (j1 = 1; j1 <= i__1; ++j1) {
55 z__[j1 + (z_dim2 << 1) * z_dim1] = 0.f;
56 i__2 = *jm;
57 for (j2 = 1; j2 <= i__2; ++j2) {
58 t1 = x[j1] - y[j2];
59 t2 = x[j1] + y[j2];
60 z__[j1 + (j2 + (z_dim2 << 1)) * z_dim1] = (sin((*jm << 1) * t1) /
61 sin(t1) - sin((*jm << 1) * t2) / sin(t2)) / (*jm << 1);
62 /* L20: */
63 }
64 z__[j1 + (*jm + (z_dim2 << 1)) * z_dim1] /= 2;
65 /* L21: */
66 }
67 i__1 = *jm;
68 for (j1 = 1; j1 <= i__1; ++j1) {
69 i__2 = *jm - 1;
70 for (j2 = 0; j2 <= i__2; ++j2) {
71 t1 = x[j1] - y[j2];
72 t2 = x[j1] + y[j2];
73 z__[j1 + (j2 + z_dim2 * 3) * z_dim1] = (sin((*jm << 1) * t1) /
74 sin(t1) + sin((*jm << 1) * t2) / sin(t2)) / (*jm << 1);
75 /* L30: */
76 }
77 z__[j1 + z_dim2 * 3 * z_dim1] /= 2;
78 /* L31: */
79 }
80 i__1 = *jm;
81 for (j1 = 1; j1 <= i__1; ++j1) {
82 z__[j1 + (z_dim2 << 2) * z_dim1] = 0.f;
83 i__2 = *jm - 1;
84 for (j2 = 1; j2 <= i__2; ++j2) {
85 t1 = x[j1] - y[j2];
86 t2 = x[j1] + y[j2];
87 z__[j1 + (j2 + (z_dim2 << 2)) * z_dim1] = (sin(((*jm << 1) + 1) *
88 t1) / sin(t1) - sin(((*jm << 1) + 1) * t2) / sin(t2)) / (*
89 jm << 1);
90 /* L40: */
91 }
92 /* L41: */
93 }
94 return 0;
95 } /* shiniz_ */
96
97