1 /*
2  * Copyright (C) 1998, 2000-2007, 2010, 2011, 2012, 2013 SINTEF ICT,
3  * Applied Mathematics, Norway.
4  *
5  * Contact information: E-mail: tor.dokken@sintef.no
6  * SINTEF ICT, Department of Applied Mathematics,
7  * P.O. Box 124 Blindern,
8  * 0314 Oslo, Norway.
9  *
10  * This file is part of SISL.
11  *
12  * SISL is free software: you can redistribute it and/or modify
13  * it under the terms of the GNU Affero General Public License as
14  * published by the Free Software Foundation, either version 3 of the
15  * License, or (at your option) any later version.
16  *
17  * SISL is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
20  * GNU Affero General Public License for more details.
21  *
22  * You should have received a copy of the GNU Affero General Public
23  * License along with SISL. If not, see
24  * <http://www.gnu.org/licenses/>.
25  *
26  * In accordance with Section 7(b) of the GNU Affero General Public
27  * License, a covered work must retain the producer line in every data
28  * file that is created or manipulated using SISL.
29  *
30  * Other Usage
31  * You can be released from the requirements of the license by purchasing
32  * a commercial license. Buying such a license is mandatory as soon as you
33  * develop commercial activities involving the SISL library without
34  * disclosing the source code of your own applications.
35  *
36  * This file may be used in accordance with the terms contained in a
37  * written agreement between you and SINTEF ICT.
38  */
39 
40 #include "sisl-copyright.h"
41 
42 /*
43  *
44  * $Id:
45  *
46  */
47 #define S1540
48 
49 
50 #include "sislP.h"
51 
52 #if defined(SISLNEEDPROTOTYPES)
53 void
s1540(double et[],int ik,int in,double ax[],int im,int ider,double ebder[],int ileft[],int * jstat)54 s1540  ( double   et[],
55 	 int      ik,
56 	 int      in,
57 	 double   ax[],
58 	 int      im,
59 	 int      ider,
60 	 double   ebder[],
61 	 int      ileft[],
62 	 int     *jstat )
63 #else
64 void
65 s1540 ( et, ik, in, ax, im, ider, ebder, ileft, jstat )
66   double   et[];
67   int      ik;
68   int      in;
69   double   ax[];
70   int      im;
71   int      ider;
72   double   ebder[];
73   int      ileft[];
74   int     *jstat;
75 #endif
76 /*
77 *********************************************************************
78 *
79 * PURPOSE:	Given a (polynomial) spline basis
80 *		(with k-regular knot vector)
81 *		and an array of parameters ax[0, ..., im-1]
82 *		calculate the basis values (and first derivatives if ider = 1)
83 *		and fill them in the array ebder[].
84 *
85 * INPUT:	ik     -  spline order
86 *		in     -  number of spline functions (only for checking ax).
87 *		ax     -  array of parameters
88 *		im     -  length of the array ax[]
89 *
90 * OUTPUT:	ebder  -  array containing the basis values
91 *		          ider = 0:  B(ax[0   ],i),...,B(ax[0   ],i+ik-1)
92 *		                     B(ax[1   ],i),...,B(ax[1   ],i+ik-1)
93 *		                     B(ax[im-1],i),...,B(ax[im-1],i+ik-1)
94 *
95 *		          ider = 1:  B (ax[0   ],i),...,B (ax[0   ],i+ik-1)
96 *		                     B'(ax[0   ],i),...,B'(ax[0   ],i+ik-1)
97 *		                     B (ax[1   ],i),...,B (ax[1   ],i+ik-1)
98 *		                     B'(ax[1   ],i),...,B'(ax[1   ],i+ik-1)
99 *						...............
100 *		                     B (ax[im-1],i),...,B (ax[im-1],i+ik-1)
101 *		                     B'(ax[im-1],i),...,B'(ax[im-1],i+ik-1)
102 *
103 *		Ie. if we have positions and derivatives, then
104 *		for each ax[i] first the positions come and then the
105 *		derivatives.
106 *              jstat    - Status messages
107 *
108 *                         = 0 : Ok.
109 *                         < 0 : Error.
110 *
111 * IMPORTANT:	ebder is filled differently by the corresponding SISL routine
112 *		(s1504):
113 *			 B (ax[0   ],i     ), B'(ax[0   ],i     ),
114 *			 B (ax[0   ],i+1   ), B'(ax[0   ],i+1   ),
115 *			        ...............
116 *		         B (ax[0   ],i+ik-1), B'(ax[0   ],i+ik-1),
117 *			 B (ax[1   ],i     ), B'(ax[1   ],i     ),
118 *			        ...............
119 *			 B (ax[im-1],i+ik-1), B'(ax[im-1],i+ik-1),
120 *
121 *
122 * WRITTEN BY: 	Geir Westgaard, SINTEF, Oslo, November 1999
123 *
124 *********************************************************************
125 */
126 {
127   int kstat = 0;      /* Local status variable.                          */
128   int kpos  = 0;      /* The position of error.                          */
129   int i, k;           /* Control variables in for loops and for stepping
130                          through arrays.                                 */
131   int size;           /* (ider+1) * ik.                                  */
132   double *eder = SISL_NULL;/* B-spline evaluationas at a single value.        */
133   double tmpeder[10]; /* meaning: tmpeder[(ider+1)*ik]
134 			 and assuming: ider <= 1 and ik <= 5             */
135 
136 
137 
138    /* Check the input. */
139 
140     if ( ider < 0 || ider > 1        ) goto err10;
141     if ( ik < 2 || ik > 5            ) goto err10;
142     if ( im  < 0                     ) goto err10;
143 
144 
145    /* Set local variables. */
146 
147    size = (ider + 1)*ik;
148 
149    eder = ebder;
150 
151 
152    if ( ider == 0 )
153    {
154       for( k = 0; k < im; k++, eder += ik )
155       {
156 	 s1220( et, ik, in, ileft + k, ax[k], ider, eder, &kstat );
157 
158 	 if ( kstat < 0 ) goto error;
159       }
160    }
161    else
162    {
163 
164       for( k = 0; k < im; k++, eder += size )
165       {
166 	 s1220( et, ik, in, ileft + k, ax[k], ider, tmpeder, &kstat );
167 
168 	 if ( kstat < 0 ) goto error;
169 
170 	 for ( i = 0; i < ik; ++i )
171 	 {
172 	    eder[i     ] = tmpeder[2*i    ];
173 	    eder[i+ik  ] = tmpeder[2*i + 1];
174 	 }
175       }
176    }
177 
178 
179 
180 
181   /* Successful computations.  */
182 
183    *jstat = 0;
184    goto out;
185 
186 /* Error in input. */
187 err10: *jstat = -10;
188 s6err( "s1540", *jstat, kpos );
189 goto out;
190 
191 /* Error in lower level routine.  */
192 error:  *jstat = kstat;
193 s6err( "s1540", *jstat, kpos );
194 goto out;
195 
196 
197 
198 
199 out: return;
200 
201 }
202