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: s9surmarch.c,v 1.2 2001-03-19 15:59:03 afr Exp $
45 *
46 */
47
48
49 #define S9SURMARCH
50
51 #include "sislP.h"
52
53 #if defined(SISLNEEDPROTOTYPES)
54 void
s9surmarch(SISLSurf * ps1,SISLSurf * ps2,double epar[],int ndir[],int ipoint,double * gpar[],int * mpar[],int * jpoint,int * jstat)55 s9surmarch(SISLSurf *ps1,SISLSurf *ps2,double epar[],int ndir[],int ipoint,
56 double *gpar[],int *mpar[],int *jpoint,int *jstat)
57 #else
58 void s9surmarch(ps1,ps2,epar,ndir,ipoint,gpar,mpar,jpoint,jstat)
59 SISLSurf *ps1;
60 SISLSurf *ps2;
61 double epar[];
62 int ndir[];
63 int ipoint;
64 double *gpar[];
65 int *mpar[];
66 int *jpoint;
67 int *jstat;
68 #endif
69 /*
70 *********************************************************************
71 *
72 *********************************************************************
73 *
74 * PURPOSE : To check which of the points in the epar array can be
75 * connected by marching. If the connection should be
76 * to internal point these points are returned. All points
77 * are assumed to lie on the boundary or close to the boundary
78 * of the patch (DEQUAL says that they lie on the boundary).
79 *
80 * When ipoint = 2, we connect always if marching succeeds!!!!!!!
81 * without testing for equality between marching results and
82 * given point no 2.
83 *
84 * INPUT : ps1 - The first surface in intersection.
85 * ps2 - The second surface in intersection.
86 * epar[2,*] - Parameter values for the intersection points along
87 * the boundary.
88 * ndir[*] - Description of the points in epar.
89 * -1 : tangent pointing out
90 * 0 : tangent parallel to boundary
91 * 1 : tangent pointing in
92 * 2 : singular point
93 * The algorithm assumes that first all points
94 * with status -1 and 1 lies first in epar.
95 * ipoint - Number of parameter values
96 *
97 *
98 * OUTPUT : gpar - Pointer to parameter values produced tor singular
99 * points. Is allocated inside the function, must
100 * be released by the calling function.
101 * mpar - Array describing how the input parameter pairs were
102 * connected to other in the input vector and to the
103 * new parameter pairs made in gpar. Is allocated
104 * inside the function, must be released by the
105 * calling function.
106 * -3 - point on closed curve
107 * -2 - internal point onc other curve
108 * -1 - Can not march from point
109 * POSITIVE VALUE - SISLPoint connected to point in epar
110 * jpoint - Number of new parameter pairs produced
111 * jstat - status messages
112 * = 2 : Only singular points, not connected
113 * = 1 : Points connected
114 * = 0 : Points not connected..
115 * < 0 : Error.
116 *
117 *
118 * METHOD :
119 *
120 *
121 * REFERENCES :
122 *
123 * WRITTEN BY : Ulf J. Krystad, SI, Oslo, Norway. August 1989
124 *
125 *********************************************************************
126 */
127 {
128 int kstat; /* Status variable */
129 int kpos=0; /* Position of error */
130 int *lpar = SISL_NULL; /* Pointer to output integer array */
131 int ki,kj;
132 int kn1,kn2,kk1,kk2; /* Surface attributes. */
133 double tstart1,tstart2,tend1,tend2; /* Surface attributes. */
134 int ksucc; /* Success indicator */
135 double tepsge=1.0; /* Not used */
136 double *spar=SISL_NULL; /* Pointer to output real array */
137 double scand1[4]; /* Result of iteration process */
138 double scand2[4]; /* Result of iteration process */
139 double *sp,*sq; /* Pointer used in loop */
140 double tdum1; /* Max knot value used in DEQUAL comparing. */
141 double tdum2; /* Max knot value used in DEQUAL comparing. */
142 double tdum3; /* Max knot value used in DEQUAL comparing. */
143 double tdum4; /* Max knot value used in DEQUAL comparing. */
144
145 /* Init */
146 kn1 = ps1->in1;
147 kn2 = ps1->in2;
148 kk1 = ps1->ik1;
149 kk2 = ps1->ik2;
150
151 tstart1 = ps1->et1[kk1-1];
152 tend1 = ps1->et1[kn1];
153 tstart2 = ps1->et2[kk2-1];
154 tend2 = ps1->et2[kn2];
155
156 tdum1 = (double)2.0*max(fabs(tstart1),fabs(tend1));
157 tdum2 = (double)2.0*max(fabs(tstart2),fabs(tend2));
158
159 kn1 = ps2->in1;
160 kn2 = ps2->in2;
161 kk1 = ps2->ik1;
162 kk2 = ps2->ik2;
163
164 tstart1 = ps2->et1[kk1-1];
165 tend1 = ps2->et1[kn1];
166 tstart2 = ps2->et2[kk2-1];
167 tend2 = ps2->et2[kn2];
168
169
170 tdum3 = (double)2.0*max(fabs(tstart1),fabs(tend1));
171 tdum4 = (double)2.0*max(fabs(tstart2),fabs(tend2));
172
173
174 /* Allocate output arrays */
175
176 if ((*mpar=newarray(2*ipoint,INT )) == SISL_NULL) goto err101;
177 if ((*gpar=newarray(8*ipoint,DOUBLE)) == SISL_NULL) goto err101;
178
179 lpar = *mpar;
180 spar = *gpar;
181
182 memcopy(spar,epar,4*ipoint,DOUBLE);
183 *jpoint = ipoint;
184
185 /* Initiate output integer array to point to no points */
186
187 for (ki=0 ; ki< 2*ipoint ; ki++) *(lpar+ki) = 0;
188
189
190 /* Loop for all input points. */
191 for (ki=0, sp=spar ; ki< ipoint-1 ; ki++, sp+=4)
192 {
193 /* Start marching from point ki */
194
195 /* Exclude points already connected and parallell points. */
196 if (lpar[ki] != 0 || ndir[ki] == 0) continue;
197
198 /* SISLPoint not marched to */
199
200 s1788(ps1,ps2,tepsge,sp,scand1,scand2,&kstat);
201 if (kstat<0) goto error;
202 if (kstat==0) goto war00;;
203
204 /* Run through remaining points to find if scand2 matches any
205 of them. If we've got only two points, we connect them.*/
206
207 ksucc = 0;
208
209 for (kj=ki+1,sq=spar+4*ki+4 ; kj<ipoint ; kj++,sq+=4)
210 {
211
212 /* SISLPoint found */
213
214 if (DEQUAL(sq[0]+tdum1,scand2[0]+tdum1) &&
215 DEQUAL(sq[1]+tdum2,scand2[1]+tdum2) &&
216 DEQUAL(sq[2]+tdum3,scand2[2]+tdum3) &&
217 DEQUAL(sq[3]+tdum4,scand2[3]+tdum4))
218
219 {
220 /* Accepted end point found */
221
222 lpar[ki] = kj+1;
223 lpar[kj] = ki+1;
224 ksucc = 1;
225 break;
226 }
227 }
228 /* If ksucc==0 then one of the searches was not successful */
229
230 if (ksucc==0) goto war00;
231
232 }
233
234 goto success;
235
236 success:
237 *jstat = 1;
238 goto out;
239
240 /* No success */
241 war00:
242 *jstat=0;
243 /* If we got only singular points, set status.*/
244 if (ndir[0] == 2) *jstat = 2;
245 goto out;
246
247 /* Error in space allocation */
248 err101:
249 *jstat = -101;
250 s6err("s9surmarch",*jstat,kpos);
251 goto out;
252
253 /* Error in lower level function */
254 error:
255 *jstat = kstat;
256 s6err("s9surmarch",*jstat,kpos);
257 goto out;
258
259 out:;
260 }
261