0023952: Improving thread-safety of intersections, approximations and other modeling...
[occt.git] / src / AdvApp2Var / AdvApp2Var_ApproxF2var.cxx
CommitLineData
b311480e 1// Copyright (c) 1999-2012 OPEN CASCADE SAS
7fd59977 2//
b311480e 3// The content of this file is subject to the Open CASCADE Technology Public
4// License Version 6.5 (the "License"). You may not use the content of this file
5// except in compliance with the License. Please obtain a copy of the License
6// at http://www.opencascade.org and read it completely before using this file.
7//
8// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
9// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
7fd59977 10//
b311480e 11// The Original Code and all software distributed under the License is
12// distributed on an "AS IS" basis, without warranty of any kind, and the
13// Initial Developer hereby disclaims all such warranties, including without
14// limitation, any warranties of merchantability, fitness for a particular
15// purpose or non-infringement. Please see the License for the specific terms
16// and conditions governing the rights and limitations under the License.
17
18// AdvApp2Var_ApproxF2var.cxx
7fd59977 19#include <math.h>
20#include <AdvApp2Var_SysBase.hxx>
21#include <AdvApp2Var_MathBase.hxx>
22#include <AdvApp2Var_Data_f2c.hxx>
23#include <AdvApp2Var_Data.hxx>
24#include <AdvApp2Var_ApproxF2var.hxx>
25
26
27static
28int mmjacpt_(const integer *ndimen,
29 const integer *ncoefu,
30 const integer *ncoefv,
31 const integer *iordru,
32 const integer *iordrv,
33 const doublereal *ptclgd,
34 doublereal *ptcaux,
35 doublereal *ptccan);
36
37
38
39static
40int mma2ce2_(integer *numdec,
41 integer *ndimen,
42 integer *nbsesp,
43 integer *ndimse,
44 integer *ndminu,
45 integer *ndminv,
46 integer *ndguli,
47 integer *ndgvli,
48 integer *ndjacu,
49 integer *ndjacv,
50 integer *iordru,
51 integer *iordrv,
52 integer *nbpntu,
53 integer *nbpntv,
54 doublereal *epsapr,
55 doublereal *sosotb,
56 doublereal *disotb,
57 doublereal *soditb,
58 doublereal *diditb,
59 doublereal *gssutb,
60 doublereal *gssvtb,
61 doublereal *xmaxju,
62 doublereal *xmaxjv,
63 doublereal *vecerr,
64 doublereal *chpair,
65 doublereal *chimpr,
66 doublereal *patjac,
67 doublereal *errmax,
68 doublereal *errmoy,
69 integer *ndegpu,
70 integer *ndegpv,
71 integer *itydec,
72 integer *iercod);
73
74static
75int mma2cfu_(integer *ndujac,
76 integer *nbpntu,
77 integer *nbpntv,
78 doublereal *sosotb,
79 doublereal *disotb,
80 doublereal *soditb,
81 doublereal *diditb,
82 doublereal *gssutb,
83 doublereal *chpair,
84 doublereal *chimpr);
85
86static
87int mma2cfv_(integer *ndvjac,
88 integer *mindgu,
89 integer *maxdgu,
90 integer *nbpntv,
91 doublereal *gssvtb,
92 doublereal *chpair,
93 doublereal *chimpr,
94 doublereal *patjac);
95
96static
97int mma2er1_(integer *ndjacu,
98 integer *ndjacv,
99 integer *ndimen,
100 integer *mindgu,
101 integer *maxdgu,
102 integer *mindgv,
103 integer *maxdgv,
104 integer *iordru,
105 integer *iordrv,
106 doublereal *xmaxju,
107 doublereal *xmaxjv,
108 doublereal *patjac,
109 doublereal *vecerr,
110 doublereal *erreur);
111
112static
113int mma2er2_(integer *ndjacu,
114 integer *ndjacv,
115 integer *ndimen,
116 integer *mindgu,
117 integer *maxdgu,
118 integer *mindgv,
119 integer *maxdgv,
120 integer *iordru,
121 integer *iordrv,
122 doublereal *xmaxju,
123 doublereal *xmaxjv,
124 doublereal *patjac,
125 doublereal *epmscut,
126 doublereal *vecerr,
127 doublereal *erreur,
128 integer *newdgu,
129 integer *newdgv);
130
131static
132int mma2moy_(integer *ndgumx,
133 integer *ndgvmx,
134 integer *ndimen,
135 integer *mindgu,
136 integer *maxdgu,
137 integer *mindgv,
138 integer *maxdgv,
139 integer *iordru,
140 integer *iordrv,
141 doublereal *patjac,
142 doublereal *errmoy);
143
144static
145int mma2ds2_(integer *ndimen,
146 doublereal *uintfn,
147 doublereal *vintfn,
41194117 148 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 149 integer *nbpntu,
150 integer *nbpntv,
151 doublereal *urootb,
152 doublereal *vrootb,
153 integer *iiuouv,
154 doublereal *sosotb,
155 doublereal *disotb,
156 doublereal *soditb,
157 doublereal *diditb,
158 doublereal *fpntab,
159 doublereal *ttable,
160 integer *iercod);
161
162
163
164
165static
166int mma1fdi_(integer *ndimen,
167 doublereal *uvfonc,
41194117 168 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 169 integer *isofav,
170 doublereal *tconst,
171 integer *nbroot,
172 doublereal *ttable,
173 integer *iordre,
174 integer *ideriv,
175 doublereal *fpntab,
176 doublereal *somtab,
177 doublereal *diftab,
178 doublereal *contr1,
179 doublereal *contr2,
180 integer *iercod);
181
182static
183int mma1cdi_(integer *ndimen,
184 integer *nbroot,
185 doublereal *rootlg,
186 integer *iordre,
187 doublereal *contr1,
188 doublereal *contr2,
189 doublereal *somtab,
190 doublereal *diftab,
191 doublereal *fpntab,
192 doublereal *hermit,
193 integer *iercod);
194static
195int mma1jak_(integer *ndimen,
196 integer *nbroot,
197 integer *iordre,
198 integer *ndgjac,
199 doublereal *somtab,
200 doublereal *diftab,
201 doublereal *cgauss,
202 doublereal *crvjac,
203 integer *iercod);
204static
205int mma1cnt_(integer *ndimen,
206 integer *iordre,
207 doublereal *contr1,
208 doublereal *contr2,
209 doublereal *hermit,
210 integer *ndgjac,
211 doublereal *crvjac);
212
213static
214int mma1fer_(integer *ndimen,
215 integer *nbsesp,
216 integer *ndimse,
217 integer *iordre,
218 integer *ndgjac,
219 doublereal *crvjac,
220 integer *ncflim,
221 doublereal *epsapr,
222 doublereal *ycvmax,
223 doublereal *errmax,
224 doublereal *errmoy,
225 integer *ncoeff,
226 integer *iercod);
227
228static
229int mma1noc_(doublereal *dfuvin,
230 integer *ndimen,
231 integer *iordre,
232 doublereal *cntrin,
233 doublereal *duvout,
234 integer *isofav,
235 integer *ideriv,
236 doublereal *cntout);
237
238
239static
240 int mmmapcoe_(integer *ndim,
241 integer *ndgjac,
242 integer *iordre,
243 integer *nbpnts,
244 doublereal *somtab,
245 doublereal *diftab,
246 doublereal *gsstab,
247 doublereal *crvjac);
248
249static
250 int mmaperm_(integer *ncofmx,
251 integer *ndim,
252 integer *ncoeff,
253 integer *iordre,
254 doublereal *crvjac,
255 integer *ncfnew,
256 doublereal *errmoy);
257
258
259#define mmapgss_1 mmapgss_
260#define mmapgs0_1 mmapgs0_
261#define mmapgs1_1 mmapgs1_
262#define mmapgs2_1 mmapgs2_
263
264//=======================================================================
265//function : mma1cdi_
266//purpose :
267//=======================================================================
268int mma1cdi_(integer *ndimen,
269 integer *nbroot,
270 doublereal *rootlg,
271 integer *iordre,
272 doublereal *contr1,
273 doublereal *contr2,
274 doublereal *somtab,
275 doublereal *diftab,
276 doublereal *fpntab,
277 doublereal *hermit,
278 integer *iercod)
279{
1ef32e96 280 integer c__1 = 1;
7fd59977 281
282 /* System generated locals */
283 integer contr1_dim1, contr1_offset, contr2_dim1, contr2_offset,
284 somtab_dim1, somtab_offset, diftab_dim1, diftab_offset,
285 fpntab_dim1, fpntab_offset, hermit_dim1, hermit_offset, i__1,
286 i__2, i__3;
41194117 287
7fd59977 288 /* Local variables */
1ef32e96
RL
289 integer nroo2, ncfhe, nd, ii, kk;
290 integer ibb, kkm, kkp;
291 doublereal bid1, bid2, bid3;
7fd59977 292
7fd59977 293/* **********************************************************************
294*/
0d969553 295/* FUNCTION : */
7fd59977 296/* ---------- */
0d969553
Y
297/* Discretisation on the parameters of interpolation polynomes */
298/* constraints of order IORDRE. */
7fd59977 299
0d969553 300/* KEYWORDS : */
7fd59977 301/* ----------- */
0d969553 302/* ALL, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
7fd59977 303
0d969553 304/* INPUT ARGUMENTS : */
7fd59977 305/* ------------------ */
0d969553
Y
306/* NDIMEN: Space dimension. */
307/* NBROOT: Number of INTERNAL discretisation parameters. */
308/* It is also the root number Legendre polynome where */
309/* the discretization is performed. */
310/* ROOTLG: Table of discretization parameters ON (-1,1). */
311/* IORDRE: Order of constraint imposed to the extremities of the iso. */
312/* = 0, the extremities of the iso are calculated */
313/* = 1, additionally, the 1st derivative in the direction */
314/* of the iso is calculated. */
315/* = 2, additionally, the 2nd derivative in the direction */
316/* of the iso is calculated. */
317/* CONTR1: Contains, if IORDRE>=0, values IORDRE+1 in TTABLE(0)
318*/
319/* (1st extremity) of derivatives of F(Uc,Ve) or F(Ue,Vc), */
320/* see below. */
321/* CONTR2: Contains, if IORDRE>=0, values IORDRE+1 in */
322/* TTABLE(NBROOT+1) (2nd extremity) of: */
323/* If ISOFAV=1, derived of order IDERIV by U, derived */
324/* ordre 0 to IORDRE by V of F(Uc,Ve) or Uc=TCONST */
325/* (fixed iso value) and Ve is the fixed extremity. */
326/* If ISOFAV=2, derivative of order IDERIV by V, derivative */
327/* of order 0 to IORDRE by U of F(Ue,Vc) or Vc=TCONST */
328/* (fixed iso value) and Ue is the fixed extremity. */
329
330/* SOMTAB: Table of NBROOT/2 sums of 2 index points */
331/* NBROOT-II+1 and II, for II = 1, NBROOT/2. */
332/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
333/* NBROOT-II+1 and II, for II = 1, NBROOT/2. */
334
335/* OUTPUT ARGUMENTS : */
7fd59977 336/* ------------------- */
0d969553
Y
337/* SOMTAB: Table of NBROOT/2 sums of 2 index points */
338/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
339/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
340/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
341/* FPNTAB: Auxiliary table. */
342/* HERMIT: Table of coeff. 2*(IORDRE+1) Hermite polynoms */
343/* of degree 2*IORDRE+1. */
344/* IERCOD: Error code, */
345/* = 0, Everythig is OK */
346/* = 1, The value of IORDRE is out of (0,2) */
347/* COMMON USED : */
7fd59977 348/* ---------------- */
349
0d969553 350/* REFERENCES CALLED : */
7fd59977 351/* ----------------------- */
352
0d969553 353/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 354/* ----------------------------------- */
0d969553
Y
355/* The results of discretization are arranged in 2 tables */
356/* SOMTAB and DIFTAB to earn time during the */
357/* calculation of coefficients of the approximation curve. */
7fd59977 358
0d969553
Y
359/* If NBROOT is uneven in SOMTAB(0,*) and DIFTAB(0,*) one stores */
360/* the values of the median root of Legendre (0.D0 in (-1,1)). */
7fd59977 361
7fd59977 362/* **********************************************************************
363*/
364
0d969553 365/* Name of the routine */
7fd59977 366
367
368 /* Parameter adjustments */
369 diftab_dim1 = *nbroot / 2 + 1;
370 diftab_offset = diftab_dim1;
371 diftab -= diftab_offset;
372 somtab_dim1 = *nbroot / 2 + 1;
373 somtab_offset = somtab_dim1;
374 somtab -= somtab_offset;
375 --rootlg;
376 hermit_dim1 = (*iordre << 1) + 2;
377 hermit_offset = hermit_dim1;
378 hermit -= hermit_offset;
379 fpntab_dim1 = *nbroot;
380 fpntab_offset = fpntab_dim1 + 1;
381 fpntab -= fpntab_offset;
382 contr2_dim1 = *ndimen;
383 contr2_offset = contr2_dim1 + 1;
384 contr2 -= contr2_offset;
385 contr1_dim1 = *ndimen;
386 contr1_offset = contr1_dim1 + 1;
387 contr1 -= contr1_offset;
388
389 /* Function Body */
390 ibb = AdvApp2Var_SysBase::mnfndeb_();
391 if (ibb >= 3) {
392 AdvApp2Var_SysBase::mgenmsg_("MMA1CDI", 7L);
393 }
394 *iercod = 0;
395
0d969553 396/* --- Recuperate 2*(IORDRE+1) coeff of 2*(IORDRE+1) of Hermite polynom ---
7fd59977 397*/
398
399 AdvApp2Var_ApproxF2var::mma1her_(iordre, &hermit[hermit_offset], iercod);
400 if (*iercod > 0) {
401 goto L9100;
402 }
403
0d969553 404/* ------------------- Discretization of Hermite polynoms -----------
7fd59977 405*/
406
407 ncfhe = (*iordre + 1) << 1;
408 i__1 = ncfhe;
409 for (ii = 1; ii <= i__1; ++ii) {
410 i__2 = *nbroot;
411 for (kk = 1; kk <= i__2; ++kk) {
412 AdvApp2Var_MathBase::mmmpocur_(&ncfhe, &c__1, &ncfhe, &hermit[ii * hermit_dim1], &
413 rootlg[kk], &fpntab[kk + ii * fpntab_dim1]);
414/* L200: */
415 }
416/* L100: */
417 }
418
0d969553 419/* ---- Discretizations of boundary polynoms are taken ----
7fd59977 420*/
421
422 nroo2 = *nbroot / 2;
423 i__1 = *ndimen;
424 for (nd = 1; nd <= i__1; ++nd) {
425 i__2 = *iordre + 1;
426 for (ii = 1; ii <= i__2; ++ii) {
427 bid1 = contr1[nd + ii * contr1_dim1];
428 bid2 = contr2[nd + ii * contr2_dim1];
429 i__3 = nroo2;
430 for (kk = 1; kk <= i__3; ++kk) {
431 kkm = nroo2 - kk + 1;
432 bid3 = bid1 * fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1] +
433 bid2 * fpntab[kkm + (ii << 1) * fpntab_dim1];
434 somtab[kk + nd * somtab_dim1] -= bid3;
435 diftab[kk + nd * diftab_dim1] += bid3;
436/* L500: */
437 }
438 i__3 = nroo2;
439 for (kk = 1; kk <= i__3; ++kk) {
440 kkp = (*nbroot + 1) / 2 + kk;
441 bid3 = bid1 * fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] +
442 bid2 * fpntab[kkp + (ii << 1) * fpntab_dim1];
443 somtab[kk + nd * somtab_dim1] -= bid3;
444 diftab[kk + nd * diftab_dim1] -= bid3;
445/* L600: */
446 }
447/* L400: */
448 }
449/* L300: */
450 }
451
0d969553 452/* ------------ Cas when discretization is done on the roots of a -----------
7fd59977 453*/
0d969553 454/* ---------- Legendre polynom of uneven degree, 0 is root --------
7fd59977 455*/
456
457 if (*nbroot % 2 == 1) {
458 i__1 = *ndimen;
459 for (nd = 1; nd <= i__1; ++nd) {
460 i__2 = *iordre + 1;
461 for (ii = 1; ii <= i__2; ++ii) {
462 bid3 = fpntab[nroo2 + 1 + ((ii << 1) - 1) * fpntab_dim1] *
463 contr1[nd + ii * contr1_dim1] + fpntab[nroo2 + 1 + (
464 ii << 1) * fpntab_dim1] * contr2[nd + ii *
465 contr2_dim1];
466/* L800: */
467 }
468 somtab[nd * somtab_dim1] -= bid3;
469 diftab[nd * diftab_dim1] -= bid3;
470/* L700: */
471 }
472 }
473
474 goto L9999;
475
476/* ------------------------------ The End -------------------------------
477*/
0d969553 478/* --> IORDRE is not in the authorized zone. */
7fd59977 479L9100:
480 *iercod = 1;
481 goto L9999;
482
483L9999:
484 if (ibb >= 3) {
485 AdvApp2Var_SysBase::mgsomsg_("MMA1CDI", 7L);
486 }
487 return 0;
488} /* mma1cdi_ */
489
490//=======================================================================
491//function : mma1cnt_
492//purpose :
493//=======================================================================
494int mma1cnt_(integer *ndimen,
495 integer *iordre,
496 doublereal *contr1,
497 doublereal *contr2,
498 doublereal *hermit,
499 integer *ndgjac,
500 doublereal *crvjac)
501{
502 /* System generated locals */
503 integer contr1_dim1, contr1_offset, contr2_dim1, contr2_offset,
504 hermit_dim1, hermit_offset, crvjac_dim1, crvjac_offset, i__1,
505 i__2, i__3;
41194117 506
7fd59977 507 /* Local variables */
1ef32e96
RL
508 integer nd, ii, jj, ibb;
509 doublereal bid;
41194117
K
510
511
7fd59977 512 /* ***********************************************************************
513 */
514
0d969553 515 /* FUNCTION : */
7fd59977 516 /* ---------- */
0d969553 517 /* Add constraint to polynom. */
7fd59977 518
519 /* MOTS CLES : */
520 /* ----------- */
0d969553 521 /* ALL,AB_SPECIFI::COURE&,APPROXIMATION,ADDITION,&CONSTRAINT */
7fd59977 522
0d969553 523 /* INPUT ARGUMENTS : */
7fd59977 524 /* -------------------- */
0d969553
Y
525 /* NDIMEN: Dimension of the space */
526 /* IORDRE: Order of constraint. */
527 /* CONTR1: pt of constraint in -1, from order 0 to IORDRE. */
528 /* CONTR2: Pt of constraint in +1, from order 0 to IORDRE. */
529 /* HERMIT: Table of Hermit polynoms of order IORDRE. */
530 /* CRVJAV: Curve of approximation in Jacobi base. */
7fd59977 531
0d969553 532 /* OUTPUT ARGUMENTS : */
7fd59977 533 /* --------------------- */
0d969553
Y
534 /* CRVJAV: Curve of approximation in Jacobi base */
535 /* to which the polynom of interpolation of constraints is added. */
7fd59977 536
0d969553 537 /* COMMON USED : */
7fd59977 538 /* ------------------ */
539
540
0d969553 541 /* REFERENCES CALLED : */
7fd59977 542 /* --------------------- */
543
544
0d969553 545/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 546/* ----------------------------------- */
547
7fd59977 548/* > */
549/* ***********************************************************************
550 */
551/* DECLARATIONS */
552/* ***********************************************************************
553 */
0d969553 554/* Name of the routine */
7fd59977 555
556/* ***********************************************************************
557 */
558/* INITIALISATIONS */
559/* ***********************************************************************
560 */
561
562 /* Parameter adjustments */
563 hermit_dim1 = (*iordre << 1) + 2;
564 hermit_offset = hermit_dim1;
565 hermit -= hermit_offset;
566 contr2_dim1 = *ndimen;
567 contr2_offset = contr2_dim1 + 1;
568 contr2 -= contr2_offset;
569 contr1_dim1 = *ndimen;
570 contr1_offset = contr1_dim1 + 1;
571 contr1 -= contr1_offset;
572 crvjac_dim1 = *ndgjac + 1;
573 crvjac_offset = crvjac_dim1;
574 crvjac -= crvjac_offset;
575
576 /* Function Body */
577 ibb = AdvApp2Var_SysBase::mnfndeb_();
578 if (ibb >= 3) {
579 AdvApp2Var_SysBase::mgenmsg_("MMA1CNT", 7L);
580 }
581
582/* ***********************************************************************
583 */
0d969553 584/* Processing */
7fd59977 585/* ***********************************************************************
586 */
587
588 i__1 = *ndimen;
589 for (nd = 1; nd <= i__1; ++nd) {
590 i__2 = (*iordre << 1) + 1;
591 for (ii = 0; ii <= i__2; ++ii) {
592 bid = 0.;
593 i__3 = *iordre + 1;
594 for (jj = 1; jj <= i__3; ++jj) {
595 bid = bid + contr1[nd + jj * contr1_dim1] *
596 hermit[ii + ((jj << 1) - 1) * hermit_dim1] +
597 contr2[nd + jj * contr2_dim1] * hermit[ii + (jj << 1) * hermit_dim1];
598 /* L300: */
599 }
600 crvjac[ii + nd * crvjac_dim1] = bid;
601 /* L200: */
602 }
603 /* L100: */
604 }
605
606/* ***********************************************************************
607 */
0d969553 608/* RETURN CALLING PROGRAM */
7fd59977 609/* ***********************************************************************
610 */
611
612 if (ibb >= 3) {
613 AdvApp2Var_SysBase::mgsomsg_("MMA1CNT", 7L);
614 }
615
616 return 0 ;
617} /* mma1cnt_ */
618
619//=======================================================================
620//function : mma1fdi_
621//purpose :
622//=======================================================================
623int mma1fdi_(integer *ndimen,
624 doublereal *uvfonc,
41194117 625 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 626 integer *isofav,
627 doublereal *tconst,
628 integer *nbroot,
629 doublereal *ttable,
630 integer *iordre,
631 integer *ideriv,
632 doublereal *fpntab,
633 doublereal *somtab,
634 doublereal *diftab,
635 doublereal *contr1,
636 doublereal *contr2,
637 integer *iercod)
638{
639 /* System generated locals */
640 integer fpntab_dim1, somtab_dim1, somtab_offset, diftab_dim1,
641 diftab_offset, contr1_dim1, contr1_offset, contr2_dim1,
642 contr2_offset, i__1, i__2;
643 doublereal d__1;
41194117 644
7fd59977 645 /* Local variables */
1ef32e96
RL
646 integer ideb, ifin, nroo2, ideru, iderv;
647 doublereal renor;
648 integer ii, nd, ibb, iim, nbp, iip;
649 doublereal bid1, bid2;
41194117 650
7fd59977 651/* **********************************************************************
652*/
653
0d969553 654/* FUNCTION : */
7fd59977 655/* ---------- */
0d969553
Y
656/* DiscretiZation of a non-polynomial function F(U,V) or of */
657/* its derivative with fixed isoparameter. */
7fd59977 658
0d969553 659/* KEYWORDS : */
7fd59977 660/* ----------- */
0d969553 661/* ALL, AB_SPECIFI::FONCTION&, DISCRETISATION, &POINT */
7fd59977 662
0d969553 663/* INPUT ARGUMENTS : */
7fd59977 664/* ------------------ */
0d969553
Y
665/* NDIMEN: Space dimension. */
666/* UVFONC: Limits of the path of definition by U and by V of the approximated function */
667/* FONCNP: The NAME of the non-polynomial function to be approximated */
668/* (external program). */
669/* ISOFAV: Fixed isoparameter for the discretization; */
670/* = 1, discretization with fixed U and variable V. */
671/* = 2, discretization with fixed V and variable U. */
672/* TCONST: Iso value is also fixed. */
673/* NBROOT: Number of INTERNAL discretization parameters. */
674/* (if there are constraints, 2 extremities should be added).
675*/
676/* This is also the root number of the Legendre polynom where */
677/* the discretization is done. */
678/* TTABLE: Table of discretization parameters and of 2 extremities */
679/* (Respectively (-1, NBROOT Legendre roots,1) */
680/* reframed within the adequate interval. */
681/* IORDRE: Order of constraint imposed on the extremities of the iso. */
682/* (If Iso-U, it is necessary to calculate the derivatives by V and vice */
7fd59977 683/* versa). */
0d969553
Y
684/* = 0, the extremities of the iso are calculated. */
685/* = 1, additionally the 1st derivative in the direction of the iso is calculated */
686/* = 2, additionally the 2nd derivative in the direction of the iso is calculated */
687/* IDERIV: Order of derivative transversal to fixed iso (If Iso-U=Uc */
688/* is fixed, the derivative of order IDERIV is discretized by U of */
689/* F(Uc,v). Same if iso-V is fixed). */
690/* Varies from 0 (positioning) to 2 (2nd derivative). */
691
692/* OUTPUT ARGUMENTS : */
7fd59977 693/* ------------------- */
0d969553
Y
694/* FPNTAB: Auxiliary table.
695 SOMTAB: Table of NBROOT/2 sums of 2 index points */
696/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
697/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
698/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
699/* CONTR1: Contains, if IORDRE>=0, values IORDRE+1 in TTABLE(0)
700*/
701/* (1st extremity) of derivatives of F(Uc,Ve) or F(Ue,Vc), */
702/* see below. */
703/* CONTR2: Contains, if IORDRE>=0, values IORDRE+1 in */
704/* TTABLE(NBROOT+1) (2nd extremity) of: */
705/* If ISOFAV=1, derived of order IDERIV by U, derived */
706/* ordre 0 to IORDRE by V of F(Uc,Ve) or Uc=TCONST */
707/* (fixed iso value) and Ve is the fixed extremity. */
708/* If ISOFAV=2, derivative of order IDERIV by V, derivative */
709/* of order 0 to IORDRE by U of F(Ue,Vc) or Vc=TCONST */
710/* (fixed iso value) and Ue is the fixed extremity. */
711/* IERCOD: Error code > 100; Pb in evaluation of FONCNP, */
712/* the returned error code is equal to error code of FONCNP + 100. */
713
714/* COMMONS USED : */
7fd59977 715/* ---------------- */
716
0d969553 717/* REFERENCES CALLED : */
7fd59977 718/* ----------------------- */
719
0d969553 720/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 721/* ----------------------------------- */
0d969553
Y
722/* The results of discretization are arranged in 2 tables */
723/* SOMTAB and DIFTAB to earn time during the */
724/* calculation of coefficients of the approximation curve. */
7fd59977 725
0d969553
Y
726/* If NBROOT is uneven in SOMTAB(0,*) and DIFTAB(0,*) one stores */
727/* the values of the median root of Legendre (0.D0 in (-1,1)). */
7fd59977 728
0d969553
Y
729/* Function F(u,v) defined in UVFONC is reparameterized in */
730/* (-1,1)x(-1,1). Then 1st and 2nd derivatives are renormalized. */
7fd59977 731
7fd59977 732/* > */
733/* **********************************************************************
734*/
735
0d969553 736/* Name of the routine */
7fd59977 737
738
739 /* Parameter adjustments */
740 uvfonc -= 3;
741 diftab_dim1 = *nbroot / 2 + 1;
742 diftab_offset = diftab_dim1;
743 diftab -= diftab_offset;
744 somtab_dim1 = *nbroot / 2 + 1;
745 somtab_offset = somtab_dim1;
746 somtab -= somtab_offset;
747 fpntab_dim1 = *ndimen;
748 --fpntab;
749 contr2_dim1 = *ndimen;
750 contr2_offset = contr2_dim1 + 1;
751 contr2 -= contr2_offset;
752 contr1_dim1 = *ndimen;
753 contr1_offset = contr1_dim1 + 1;
754 contr1 -= contr1_offset;
755
756 /* Function Body */
757 ibb = AdvApp2Var_SysBase::mnfndeb_();
758 if (ibb >= 3) {
759 AdvApp2Var_SysBase::mgenmsg_("MMA1FDI", 7L);
760 }
761 *iercod = 0;
762
0d969553 763/* --------------- Definition of the nb of points to calculate --------------
7fd59977 764*/
0d969553 765/* --> If constraints, the limits are also taken */
7fd59977 766 if (*iordre >= 0) {
767 ideb = 0;
768 ifin = *nbroot + 1;
0d969553 769/* --> Otherwise, only Legendre roots (reframed) are used
7fd59977 770. */
771 } else {
772 ideb = 1;
773 ifin = *nbroot;
774 }
0d969553 775/* --> Nb of point to calculate. */
7fd59977 776 nbp = ifin - ideb + 1;
777 nroo2 = *nbroot / 2;
778
0d969553 779/* --------------- Determination of the order of global derivation --------
7fd59977 780*/
0d969553
Y
781/* --> ISOFAV takes only values 1 or 2. */
782/* if Iso-U, derive by U of order IDERIV */
7fd59977 783 if (*isofav == 1) {
784 ideru = *ideriv;
785 iderv = 0;
786 d__1 = (uvfonc[4] - uvfonc[3]) / 2.;
787 renor = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
0d969553 788/* if Iso-V, derive by V of order IDERIV */
7fd59977 789 } else {
790 ideru = 0;
791 iderv = *ideriv;
792 d__1 = (uvfonc[6] - uvfonc[5]) / 2.;
793 renor = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
794 }
795
0d969553 796/* ----------- Discretization on roots of the ---------------
7fd59977 797*/
0d969553 798/* ---------------------- Legendre polynom of degree NBROOT -------------------
7fd59977 799*/
800
fadcea2c 801 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (ndimen,
7fd59977 802 &uvfonc[3],
803 &uvfonc[5],
804 isofav,
805 tconst,
806 &nbp,
807 &ttable[ideb],
808 &ideru,
809 &iderv,
810 &fpntab[ideb * fpntab_dim1 + 1],
811 iercod);
812 if (*iercod > 0) {
813 goto L9999;
814 }
815 i__1 = *ndimen;
816 for (nd = 1; nd <= i__1; ++nd) {
817 i__2 = nroo2;
818 for (ii = 1; ii <= i__2; ++ii) {
819 iip = (*nbroot + 1) / 2 + ii;
820 iim = nroo2 - ii + 1;
821 bid1 = fpntab[nd + iim * fpntab_dim1];
822 bid2 = fpntab[nd + iip * fpntab_dim1];
823 somtab[ii + nd * somtab_dim1] = renor * (bid2 + bid1);
824 diftab[ii + nd * diftab_dim1] = renor * (bid2 - bid1);
825/* L200: */
826 }
827/* L100: */
828 }
829
0d969553 830/* ------------ Case when discretisation is done on roots of a ----
7fd59977 831*/
0d969553 832/* ---------- Legendre polynom of uneven degree, 0 is root --------
7fd59977 833*/
834
835 if (*nbroot % 2 == 1) {
836 i__1 = *ndimen;
837 for (nd = 1; nd <= i__1; ++nd) {
838 somtab[nd * somtab_dim1] = renor * fpntab[nd + (nroo2 + 1) *
839 fpntab_dim1];
840 diftab[nd * diftab_dim1] = renor * fpntab[nd + (nroo2 + 1) *
841 fpntab_dim1];
842/* L300: */
843 }
844 } else {
845 i__1 = *ndimen;
846 for (nd = 1; nd <= i__1; ++nd) {
847 somtab[nd * somtab_dim1] = 0.;
848 diftab[nd * diftab_dim1] = 0.;
849 }
850 }
851
852
0d969553 853/* --------------------- Take into account constraints ----------------
7fd59977 854*/
855
856 if (*iordre >= 0) {
0d969553 857/* --> Recover already calculated extremities. */
7fd59977 858 i__1 = *ndimen;
859 for (nd = 1; nd <= i__1; ++nd) {
860 contr1[nd + contr1_dim1] = renor * fpntab[nd];
861 contr2[nd + contr2_dim1] = renor * fpntab[nd + (*nbroot + 1) *
862 fpntab_dim1];
863/* L400: */
864 }
0d969553 865/* --> Nb of points to calculate/call to FONCNP */
7fd59977 866 nbp = 1;
0d969553 867/* If Iso-U, derive by V till order IORDRE */
7fd59977 868 if (*isofav == 1) {
0d969553 869/* --> Factor of normalisation 1st derivative. */
7fd59977 870 bid1 = (uvfonc[6] - uvfonc[5]) / 2.;
871 i__1 = *iordre;
872 for (iderv = 1; iderv <= i__1; ++iderv) {
fadcea2c
RL
873 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
874 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
875 nbp, ttable, &ideru, &iderv, &contr1[(iderv + 1) *
7fd59977 876 contr1_dim1 + 1], iercod);
877 if (*iercod > 0) {
878 goto L9999;
879 }
880/* L500: */
881 }
882 i__1 = *iordre;
883 for (iderv = 1; iderv <= i__1; ++iderv) {
fadcea2c
RL
884 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
885 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
886 nbp, &ttable[*nbroot + 1], &ideru, &iderv, &contr2[(
7fd59977 887 iderv + 1) * contr2_dim1 + 1], iercod);
888 if (*iercod > 0) {
889 goto L9999;
890 }
891/* L510: */
892 }
0d969553 893/* If Iso-V, derive by U till order IORDRE */
7fd59977 894 } else {
0d969553 895/* --> Factor of normalization 1st derivative. */
7fd59977 896 bid1 = (uvfonc[4] - uvfonc[3]) / 2.;
897 i__1 = *iordre;
898 for (ideru = 1; ideru <= i__1; ++ideru) {
fadcea2c
RL
899 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
900 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
901 nbp, ttable, &ideru, &iderv, &contr1[(ideru + 1) *
7fd59977 902 contr1_dim1 + 1], iercod);
903 if (*iercod > 0) {
904 goto L9999;
905 }
906/* L600: */
907 }
908 i__1 = *iordre;
909 for (ideru = 1; ideru <= i__1; ++ideru) {
fadcea2c
RL
910 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
911 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
912 nbp, &ttable[*nbroot + 1], &ideru, &iderv, &contr2[(
7fd59977 913 ideru + 1) * contr2_dim1 + 1], iercod);
914 if (*iercod > 0) {
915 goto L9999;
916 }
917/* L610: */
918 }
919 }
920
0d969553 921/* ------------------------- Normalization of derivatives -------------
7fd59977 922---- */
0d969553 923/* (The function is redefined on (-1,1)*(-1,1)) */
7fd59977 924 bid2 = renor;
925 i__1 = *iordre;
926 for (ii = 1; ii <= i__1; ++ii) {
927 bid2 = bid1 * bid2;
928 i__2 = *ndimen;
929 for (nd = 1; nd <= i__2; ++nd) {
930 contr1[nd + (ii + 1) * contr1_dim1] *= bid2;
931 contr2[nd + (ii + 1) * contr2_dim1] *= bid2;
932/* L710: */
933 }
934/* L700: */
935 }
936 }
937
938/* ------------------------------ The end -------------------------------
939*/
940
941L9999:
942 if (*iercod > 0) {
943 *iercod += 100;
944 AdvApp2Var_SysBase::maermsg_("MMA1FDI", iercod, 7L);
945 }
946 if (ibb >= 3) {
947 AdvApp2Var_SysBase::mgsomsg_("MMA1FDI", 7L);
948 }
949 return 0;
950} /* mma1fdi_ */
951
952//=======================================================================
953//function : mma1fer_
954//purpose :
955//=======================================================================
956int mma1fer_(integer *,//ndimen,
957 integer *nbsesp,
958 integer *ndimse,
959 integer *iordre,
960 integer *ndgjac,
961 doublereal *crvjac,
962 integer *ncflim,
963 doublereal *epsapr,
964 doublereal *ycvmax,
965 doublereal *errmax,
966 doublereal *errmoy,
967 integer *ncoeff,
968 integer *iercod)
969{
970 /* System generated locals */
971 integer crvjac_dim1, crvjac_offset, i__1, i__2;
41194117 972
7fd59977 973 /* Local variables */
1ef32e96
RL
974 integer idim, ncfja, ncfnw, ndses, ii, kk, ibb, ier;
975 integer nbr0;
41194117
K
976
977
7fd59977 978/* ***********************************************************************
979 */
980
0d969553 981/* FUNCTION : */
7fd59977 982/* ---------- */
0d969553 983/* Calculate the degree and the errors of approximation of a border. */
7fd59977 984
0d969553 985/* KEYWORDS : */
7fd59977 986/* ----------- */
987/* TOUS,AB_SPECIFI :: COURBE&,TRONCATURE, &PRECISION */
988
0d969553 989/* INPUT ARGUMENTS : */
7fd59977 990/* -------------------- */
7fd59977 991
0d969553
Y
992/* NDIMEN: Total Dimension of the space (sum of dimensions of sub-spaces) */
993/* NBSESP: Number of "independent" sub-spaces. */
994/* NDIMSE: Table of dimensions of sub-spaces. */
995/* IORDRE: Order of constraint at the extremities of the border */
996/* -1 = no constraints, */
997/* 0 = constraints of passage to limits (i.e. C0), */
998/* 1 = C0 + constraintes of 1st derivatives (i.e. C1), */
999/* 2 = C1 + constraintes of 2nd derivatives (i.e. C2). */
1000/* NDGJAC: Degree of development in series to use for the calculation
1001/* in the base of Jacobi. */
1002/* CRVJAC: Table of coeff. of the curve of approximation in the */
1003/* base of Jacobi. */
1004/* NCFLIM: Max number of coeff of the polynomial curve */
1005/* of approximation (should be above or equal to */
1006/* 2*IORDRE+2 and below or equal to 50). */
1007/* EPSAPR: Table of errors of approximations that cannot be passed, */
1008/* sub-space by sub-space. */
1009
1010/* OUTPUT ARGUMENTS : */
7fd59977 1011/* --------------------- */
0d969553
Y
1012/* YCVMAX: Auxiliary Table. */
1013/* ERRMAX: Table of errors (sub-space by sub-space) */
1014/* MAXIMUM made in the approximation of FONCNP by */
7fd59977 1015/* COURBE. */
0d969553
Y
1016/* ERRMOY: Table of errors (sub-space by sub-space) */
1017/* AVERAGE made in the approximation of FONCNP by */
7fd59977 1018/* COURBE. */
0d969553
Y
1019/* NCOEFF: Number of significative coeffs. of the calculated "curve". */
1020/* IERCOD: Error code */
7fd59977 1021/* = 0, ok, */
0d969553
Y
1022/* =-1, warning, required tolerance can't be */
1023/* met with coefficients NFCLIM. */
1024/* = 1, order of constraints (IORDRE) is not within authorised values */
1025/*
7fd59977 1026
0d969553 1027/* COMMONS USED : */
7fd59977 1028/* ------------------ */
1029
0d969553 1030/* REFERENCES CALLED : */
7fd59977 1031/* --------------------- */
1032
0d969553 1033/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 1034/* ----------------------------------- */
7fd59977 1035/* > */
1036/* **********************************************************************
1037*/
1038
0d969553 1039/* Name of the routine */
7fd59977 1040
1041
1042 /* Parameter adjustments */
1043 --ycvmax;
1044 --errmoy;
1045 --errmax;
1046 --epsapr;
1047 --ndimse;
1048 crvjac_dim1 = *ndgjac + 1;
1049 crvjac_offset = crvjac_dim1;
1050 crvjac -= crvjac_offset;
1051
1052 /* Function Body */
1053 ibb = AdvApp2Var_SysBase::mnfndeb_();
1054 if (ibb >= 3) {
1055 AdvApp2Var_SysBase::mgenmsg_("MMA1FER", 7L);
1056 }
1057 *iercod = 0;
1058 idim = 1;
1059 *ncoeff = 0;
1060 ncfja = *ndgjac + 1;
1061
0d969553 1062/* ------------ Calculate the degree of the curve and of the Max error --------
7fd59977 1063*/
0d969553 1064/* -------------- of approximation for all sub-spaces --------
7fd59977 1065*/
1066
1067 i__1 = *nbsesp;
1068 for (ii = 1; ii <= i__1; ++ii) {
1069 ndses = ndimse[ii];
1070
0d969553 1071/* ------------ cutting of coeff. and calculation of Max error -------
7fd59977 1072---- */
1073
1074 AdvApp2Var_MathBase::mmtrpjj_(&ncfja, &ndses, &ncfja, &epsapr[ii], iordre, &crvjac[idim *
1075 crvjac_dim1], &ycvmax[1], &errmax[ii], &ncfnw);
1076
1077/* ******************************************************************
1078**** */
0d969553 1079/* ------------- If precision OK, calculate the average error -------
7fd59977 1080---- */
1081/* ******************************************************************
1082**** */
1083
1084 if (ncfnw <= *ncflim) {
1085 mmaperm_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1086 crvjac_dim1], &ncfnw, &errmoy[ii]);
41194117 1087 *ncoeff = advapp_max(ncfnw,*ncoeff);
7fd59977 1088
0d969553 1089/* ------------- Set the declined coefficients to 0.D0 -----------
7fd59977 1090-------- */
1091
1092 nbr0 = *ncflim - ncfnw;
1093 if (nbr0 > 0) {
1094 i__2 = ndses;
1095 for (kk = 1; kk <= i__2; ++kk) {
1096 AdvApp2Var_SysBase::mvriraz_(&nbr0,
fadcea2c 1097 &crvjac[ncfnw + (idim + kk - 1) * crvjac_dim1]);
7fd59977 1098/* L200: */
1099 }
1100 }
1101 } else {
1102
1103/* **************************************************************
1104******** */
0d969553 1105/* ------------------- If required precision can't be reached----
7fd59977 1106-------- */
1107/* **************************************************************
1108******** */
1109
1110 *iercod = -1;
1111
0d969553 1112/* ------------------------- calculate the Max error ------------
7fd59977 1113-------- */
1114
1115 AdvApp2Var_MathBase::mmaperx_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1116 crvjac_dim1], ncflim, &ycvmax[1], &errmax[ii], &ier);
1117 if (ier > 0) {
1118 goto L9100;
1119 }
1120
0d969553 1121/* -------------------- nb of coeff to be returned -------------
7fd59977 1122-------- */
1123
1124 *ncoeff = *ncflim;
1125
0d969553 1126/* ------------------- and calculation of the average error ----
7fd59977 1127-------- */
1128
1129 mmaperm_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1130 crvjac_dim1], ncflim, &errmoy[ii]);
1131 }
1132 idim += ndses;
1133/* L100: */
1134 }
1135
1136 goto L9999;
1137
1138/* ------------------------------ The end -------------------------------
1139*/
0d969553 1140/* --> The order of constraints is not within autorized values. */
7fd59977 1141L9100:
1142 *iercod = 1;
1143 goto L9999;
1144
1145L9999:
1146 if (*iercod != 0) {
1147 AdvApp2Var_SysBase::maermsg_("MMA1FER", iercod, 7L);
1148 }
1149 if (ibb >= 3) {
1150 AdvApp2Var_SysBase::mgsomsg_("MMA1FER", 7L);
1151 }
1152 return 0;
1153} /* mma1fer_ */
1154
1155
1156//=======================================================================
1157//function : mma1her_
1158//purpose :
1159//=======================================================================
1160int AdvApp2Var_ApproxF2var::mma1her_(const integer *iordre,
1161 doublereal *hermit,
1162 integer *iercod)
1163{
1164 /* System generated locals */
1165 integer hermit_dim1, hermit_offset;
41194117 1166
7fd59977 1167 /* Local variables */
1ef32e96 1168 integer ibb;
41194117 1169
7fd59977 1170
1171
1172/* **********************************************************************
1173*/
1174
0d969553 1175/* FUNCTION : */
7fd59977 1176/* ---------- */
0d969553
Y
1177/* Calculate 2*(IORDRE+1) Hermit polynoms of degree 2*IORDRE+1 */
1178/* on (-1,1) */
7fd59977 1179
0d969553 1180/* KEYWORDS : */
7fd59977 1181/* ----------- */
0d969553 1182/* ALL, AB_SPECIFI::CONTRAINTE&, INTERPOLATION, &POLYNOME */
7fd59977 1183
0d969553 1184/* INPUT ARGUMENTS : */
7fd59977 1185/* ------------------ */
0d969553
Y
1186/* IORDRE: Order of constraint. */
1187/* = 0, Polynom of interpolation of order C0 on (-1,1). */
1188/* = 1, Polynom of interpolation of order C0 and C1 on (-1,1). */
1189/* = 2, Polynom of interpolation of order C0, C1 and C2 on (-1,1).
7fd59977 1190*/
1191
0d969553 1192/* OUTPUT ARGUMENTS : */
7fd59977 1193/* ------------------- */
0d969553
Y
1194/* HERMIT: Table of 2*IORDRE+2 coeff. of each of 2*(IORDRE+1) */
1195/* HERMIT polynom. */
1196/* IERCOD: Error code, */
7fd59977 1197/* = 0, Ok */
0d969553
Y
1198/* = 1, required order of constraint is not managed here. */
1199/* COMMONS USED : */
7fd59977 1200/* ---------------- */
1201
0d969553 1202/* REFERENCES CALLED : */
7fd59977 1203/* ----------------------- */
1204
0d969553 1205/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 1206/* ----------------------------------- */
0d969553
Y
1207/* The part of HERMIT(*,2*i+j) table where j=1 or 2 and i=0 to IORDRE,
1208/* contains the coefficients of the polynom of degree 2*IORDRE+1 */
1209/* such as ALL values in -1 and in +1 of this polynom and its */
1210/* derivatives till order of derivation IORDRE are NULL, */
1211/* EXCEPT for the derivative of order i: */
1212/* - valued 1 in -1 if j=1 */
1213/* - valued 1 in +1 if j=2. */
7fd59977 1214/* > */
1215/* **********************************************************************
1216*/
1217
0d969553 1218/* Name of the routine */
7fd59977 1219
1220
1221 /* Parameter adjustments */
1222 hermit_dim1 = (*iordre + 1) << 1;
1223 hermit_offset = hermit_dim1 + 1;
1224 hermit -= hermit_offset;
1225
1226 /* Function Body */
1227 ibb = AdvApp2Var_SysBase::mnfndeb_();
1228 if (ibb >= 3) {
1229 AdvApp2Var_SysBase::mgenmsg_("MMA1HER", 7L);
1230 }
1231 *iercod = 0;
1232
0d969553 1233/* --- Recover (IORDRE+2) coeff of 2*(IORDRE+1) Hermit polynoms --
7fd59977 1234*/
1235
1236 if (*iordre == 0) {
1237 hermit[hermit_dim1 + 1] = .5;
1238 hermit[hermit_dim1 + 2] = -.5;
1239
1240 hermit[(hermit_dim1 << 1) + 1] = .5;
1241 hermit[(hermit_dim1 << 1) + 2] = .5;
1242 } else if (*iordre == 1) {
1243 hermit[hermit_dim1 + 1] = .5;
1244 hermit[hermit_dim1 + 2] = -.75;
1245 hermit[hermit_dim1 + 3] = 0.;
1246 hermit[hermit_dim1 + 4] = .25;
1247
1248 hermit[(hermit_dim1 << 1) + 1] = .5;
1249 hermit[(hermit_dim1 << 1) + 2] = .75;
1250 hermit[(hermit_dim1 << 1) + 3] = 0.;
1251 hermit[(hermit_dim1 << 1) + 4] = -.25;
1252
1253 hermit[hermit_dim1 * 3 + 1] = .25;
1254 hermit[hermit_dim1 * 3 + 2] = -.25;
1255 hermit[hermit_dim1 * 3 + 3] = -.25;
1256 hermit[hermit_dim1 * 3 + 4] = .25;
1257
1258 hermit[(hermit_dim1 << 2) + 1] = -.25;
1259 hermit[(hermit_dim1 << 2) + 2] = -.25;
1260 hermit[(hermit_dim1 << 2) + 3] = .25;
1261 hermit[(hermit_dim1 << 2) + 4] = .25;
1262 } else if (*iordre == 2) {
1263 hermit[hermit_dim1 + 1] = .5;
1264 hermit[hermit_dim1 + 2] = -.9375;
1265 hermit[hermit_dim1 + 3] = 0.;
1266 hermit[hermit_dim1 + 4] = .625;
1267 hermit[hermit_dim1 + 5] = 0.;
1268 hermit[hermit_dim1 + 6] = -.1875;
1269
1270 hermit[(hermit_dim1 << 1) + 1] = .5;
1271 hermit[(hermit_dim1 << 1) + 2] = .9375;
1272 hermit[(hermit_dim1 << 1) + 3] = 0.;
1273 hermit[(hermit_dim1 << 1) + 4] = -.625;
1274 hermit[(hermit_dim1 << 1) + 5] = 0.;
1275 hermit[(hermit_dim1 << 1) + 6] = .1875;
1276
1277 hermit[hermit_dim1 * 3 + 1] = .3125;
1278 hermit[hermit_dim1 * 3 + 2] = -.4375;
1279 hermit[hermit_dim1 * 3 + 3] = -.375;
1280 hermit[hermit_dim1 * 3 + 4] = .625;
1281 hermit[hermit_dim1 * 3 + 5] = .0625;
1282 hermit[hermit_dim1 * 3 + 6] = -.1875;
1283
1284 hermit[(hermit_dim1 << 2) + 1] = -.3125;
1285 hermit[(hermit_dim1 << 2) + 2] = -.4375;
1286 hermit[(hermit_dim1 << 2) + 3] = .375;
1287 hermit[(hermit_dim1 << 2) + 4] = .625;
1288 hermit[(hermit_dim1 << 2) + 5] = -.0625;
1289 hermit[(hermit_dim1 << 2) + 6] = -.1875;
1290
1291 hermit[hermit_dim1 * 5 + 1] = .0625;
1292 hermit[hermit_dim1 * 5 + 2] = -.0625;
1293 hermit[hermit_dim1 * 5 + 3] = -.125;
1294 hermit[hermit_dim1 * 5 + 4] = .125;
1295 hermit[hermit_dim1 * 5 + 5] = .0625;
1296 hermit[hermit_dim1 * 5 + 6] = -.0625;
1297
1298 hermit[hermit_dim1 * 6 + 1] = .0625;
1299 hermit[hermit_dim1 * 6 + 2] = .0625;
1300 hermit[hermit_dim1 * 6 + 3] = -.125;
1301 hermit[hermit_dim1 * 6 + 4] = -.125;
1302 hermit[hermit_dim1 * 6 + 5] = .0625;
1303 hermit[hermit_dim1 * 6 + 6] = .0625;
1304 } else {
1305 *iercod = 1;
1306 }
1307
1308/* ------------------------------ The End -------------------------------
1309*/
1310
1311 AdvApp2Var_SysBase::maermsg_("MMA1HER", iercod, 7L);
1312 if (ibb >= 3) {
1313 AdvApp2Var_SysBase::mgsomsg_("MMA1HER", 7L);
1314 }
1315 return 0;
1316} /* mma1her_ */
1317//=======================================================================
1318//function : mma1jak_
1319//purpose :
1320//=======================================================================
1321int mma1jak_(integer *ndimen,
1322 integer *nbroot,
1323 integer *iordre,
1324 integer *ndgjac,
1325 doublereal *somtab,
1326 doublereal *diftab,
1327 doublereal *cgauss,
1328 doublereal *crvjac,
1329 integer *iercod)
1330{
1331 /* System generated locals */
1332 integer somtab_dim1, somtab_offset, diftab_dim1, diftab_offset,
1333 crvjac_dim1, crvjac_offset, cgauss_dim1;
41194117 1334
7fd59977 1335 /* Local variables */
1ef32e96 1336 integer ibb;
7fd59977 1337
1338/* **********************************************************************
1339*/
1340
0d969553 1341/* FUNCTION : */
7fd59977 1342/* ---------- */
0d969553
Y
1343/* Calculate the curve of approximation of a non-polynomial function */
1344/* in the base of Jacobi. */
7fd59977 1345
0d969553 1346/* KEYWORDS : */
7fd59977 1347/* ----------- */
0d969553 1348/* FUNCTION,DISCRETISATION,APPROXIMATION,CONSTRAINT,CURVE,JACOBI */
7fd59977 1349
0d969553 1350/* INPUT ARGUMENTS : */
7fd59977 1351/* ------------------ */
0d969553
Y
1352/* NDIMEN: Total dimension of the space (sum of dimensions */
1353/* of sub-spaces) */
1354/* NBROOT: Nb of points of discretization of the iso, extremities not
1355/* included. */
1356/* IORDRE: Order of constraint at the extremities of the boundary */
1357/* -1 = no constraints, */
1358/* 0 = constraints of passage of limits (i.e. C0), */
1359/* 1 = C0 + constraints of 1st derivatives (i.e. C1), */
1360/* 2 = C1 + constraints of 2nd derivatives (i.e. C2). */
1361/* NDGJAC: Degree of development in series to be used for calculation in the
1362/* base of Jacobi. */
1363
1364/* OUTPUT ARGUMENTS : */
7fd59977 1365/* ------------------- */
0d969553
Y
1366/* CRVJAC : Curve of approximation of FONCNP with (eventually) */
1367/* taking into account of constraints at the extremities. */
1368/* This curve is of degree NDGJAC. */
1369/* IERCOD : Error code : */
1370/* 0 = All is ok. */
1371/* 33 = Pb to return data of du block data */
1372/* of coeff. of integration by GAUSS method. */
1373/* by program MMAPPTT. */
1374
1375/* COMMONS USED : */
7fd59977 1376/* ---------------- */
1377
0d969553 1378/* REFERENCES CALLED : */
7fd59977 1379/* ----------------------- */
7fd59977 1380/* > */
1381/* **********************************************************************
1382*/
1383
0d969553 1384/* Name of the routine */
7fd59977 1385
1386 /* Parameter adjustments */
1387 diftab_dim1 = *nbroot / 2 + 1;
1388 diftab_offset = diftab_dim1;
1389 diftab -= diftab_offset;
1390 somtab_dim1 = *nbroot / 2 + 1;
1391 somtab_offset = somtab_dim1;
1392 somtab -= somtab_offset;
1393 crvjac_dim1 = *ndgjac + 1;
1394 crvjac_offset = crvjac_dim1;
1395 crvjac -= crvjac_offset;
1396 cgauss_dim1 = *nbroot / 2 + 1;
1397
1398 /* Function Body */
1399 ibb = AdvApp2Var_SysBase::mnfndeb_();
1400 if (ibb >= 2) {
1401 AdvApp2Var_SysBase::mgenmsg_("MMA1JAK", 7L);
1402 }
1403 *iercod = 0;
1404
0d969553 1405/* ----------------- Recover coeffs of integration by Gauss -----------
7fd59977 1406*/
1407
1408 AdvApp2Var_ApproxF2var::mmapptt_(ndgjac, nbroot, iordre, cgauss, iercod);
1409 if (*iercod > 0) {
1410 *iercod = 33;
1411 goto L9999;
1412 }
1413
0d969553 1414/* --------------- Calculate the curve in the base of Jacobi -----------
7fd59977 1415*/
1416
1417 mmmapcoe_(ndimen, ndgjac, iordre, nbroot, &somtab[somtab_offset], &diftab[
1418 diftab_offset], cgauss, &crvjac[crvjac_offset]);
1419
1420/* ------------------------------ The End -------------------------------
1421*/
1422
1423L9999:
1424 if (*iercod != 0) {
1425 AdvApp2Var_SysBase::maermsg_("MMA1JAK", iercod, 7L);
1426 }
1427 if (ibb >= 2) {
1428 AdvApp2Var_SysBase::mgsomsg_("MMA1JAK", 7L);
1429 }
1430 return 0;
1431} /* mma1jak_ */
1432
1433//=======================================================================
1434//function : mma1noc_
1435//purpose :
1436//=======================================================================
1437int mma1noc_(doublereal *dfuvin,
1438 integer *ndimen,
1439 integer *iordre,
1440 doublereal *cntrin,
1441 doublereal *duvout,
1442 integer *isofav,
1443 integer *ideriv,
1444 doublereal *cntout)
1445{
1446 /* System generated locals */
1447 integer i__1;
1448 doublereal d__1;
41194117 1449
7fd59977 1450 /* Local variables */
1ef32e96
RL
1451 doublereal rider, riord;
1452 integer nd, ibb;
1453 doublereal bid;
7fd59977 1454/* **********************************************************************
1455*/
1456
0d969553 1457/* FUNCTION : */
7fd59977 1458/* ---------- */
0d969553
Y
1459/* Normalization of constraints of derivatives, defined on DFUVIN */
1460/* on block DUVOUT. */
7fd59977 1461
0d969553 1462/* KEYWORDS : */
7fd59977 1463/* ----------- */
0d969553 1464/* ALL, AB_SPECIFI::VECTEUR&,DERIVEE&,NORMALISATION,&VECTEUR */
7fd59977 1465
0d969553 1466/* INPUT ARGUMENTS : */
7fd59977 1467/* ------------------ */
0d969553 1468/* DFUVIN: Limits of the block of definition by U and by V where
7fd59977 1469*/
0d969553
Y
1470/* constraints CNTRIN are defined. */
1471/* NDIMEN: Dimension of the space. */
1472/* IORDRE: Order of constraint imposed at the extremities of the iso. */
1473/* (if Iso-U, it is necessary to calculate derivatives by V and vice */
7fd59977 1474/* versa). */
0d969553
Y
1475/* = 0, the extremities of the iso are calculated */
1476/* = 1, additionally the 1st derivative in the direction */
1477/* of the iso is calculated */
1478/* = 2, additionally the 2nd derivative in the direction */
1479/* of the iso is calculated */
1480/* CNTRIN: Contains, if IORDRE>=0, IORDRE+1 derivatives */
1481/* of order IORDRE of F(Uc,v) or of F(u,Vc), following the */
1482/* value of ISOFAV, RENORMALIZED by u and v in (-1,1). */
1483/* DUVOUT: Limits of the block of definition by U and by V where the */
1484/* constraints CNTOUT will be defined. */
1485/* ISOFAV: Isoparameter fixed for the discretization; */
1486/* = 1, discretization with fixed U=Uc and variable V. */
1487/* = 2, discretization with fixed V=Vc and variable U. */
7fd59977 1488/* IDERIV: Ordre de derivee transverse a l'iso fixee (Si Iso-U=Uc */
0d969553
Y
1489/* is fixed, the derivative of order IDERIV is discretized by U */
1490/* of F(Uc,v). The same if iso-V is fixed). */
1491/* Varies from (positioning) to 2 (2nd derivative). */
7fd59977 1492
0d969553 1493/* OUTPUT ARGUMENTS : */
7fd59977 1494/* ------------------- */
0d969553
Y
1495/* CNTOUT: Contains, if IORDRE>=0, IORDRE+1 derivatives */
1496/* of order IORDRE of F(Uc,v) or of F(u,Vc), depending on the */
1497/* value of ISOFAV, RENORMALIZED for u and v in DUVOUT. */
7fd59977 1498
0d969553 1499/* COMMONS USED : */
7fd59977 1500/* ---------------- */
1501
0d969553
Y
1502/* REFERENCES CALLED : */
1503/* --------------------- */
7fd59977 1504
0d969553
Y
1505/* DESCRIPTION/NOTES/LIMITATIONS : */
1506/* ------------------------------- */
1507/* CNTRIN can be an output/input argument, */
1508/* so the call: */
7fd59977 1509
1510/* CALL MMA1NOC(DFUVIN,NDIMEN,IORDRE,CNTRIN,DUVOUT */
1511/* 1 ,ISOFAV,IDERIV,CNTRIN) */
1512
0d969553 1513/* is correct. */
7fd59977 1514/* > */
1515/* **********************************************************************
1516*/
1517
0d969553 1518/* Name of the routine */
7fd59977 1519
1520
1521 /* Parameter adjustments */
1522 dfuvin -= 3;
1523 --cntout;
1524 --cntrin;
1525 duvout -= 3;
1526
1527 /* Function Body */
1528 ibb = AdvApp2Var_SysBase::mnfndeb_();
1529 if (ibb >= 3) {
1530 AdvApp2Var_SysBase::mgenmsg_("MMA1NOC", 7L);
1531 }
1532
0d969553 1533/* --------------- Determination of coefficients of normalization -------
7fd59977 1534 */
1535
1536 if (*isofav == 1) {
1537 d__1 = (dfuvin[4] - dfuvin[3]) / (duvout[4] - duvout[3]);
1538 rider = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
1539 d__1 = (dfuvin[6] - dfuvin[5]) / (duvout[6] - duvout[5]);
1540 riord = AdvApp2Var_MathBase::pow__di(&d__1, iordre);
1541
1542 } else {
1543 d__1 = (dfuvin[6] - dfuvin[5]) / (duvout[6] - duvout[5]);
1544 rider = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
1545 d__1 = (dfuvin[4] - dfuvin[3]) / (duvout[4] - duvout[3]);
1546 riord = AdvApp2Var_MathBase::pow__di(&d__1, iordre);
1547 }
1548
0d969553 1549/* ------------- Renormalization of the vector of constraint ---------------
7fd59977 1550*/
1551
1552 bid = rider * riord;
1553 i__1 = *ndimen;
1554 for (nd = 1; nd <= i__1; ++nd) {
1555 cntout[nd] = bid * cntrin[nd];
1556/* L100: */
1557 }
1558
1559/* ------------------------------ The end -------------------------------
1560*/
1561
1562 if (ibb >= 3) {
1563 AdvApp2Var_SysBase::mgsomsg_("MMA1NOC", 7L);
1564 }
1565 return 0;
1566} /* mma1noc_ */
1567
1568//=======================================================================
1569//function : mma1nop_
1570//purpose :
1571//=======================================================================
1572int mma1nop_(integer *nbroot,
1573 doublereal *rootlg,
1574 doublereal *uvfonc,
1575 integer *isofav,
1576 doublereal *ttable,
1577 integer *iercod)
1578
1579{
1580 /* System generated locals */
1581 integer i__1;
41194117 1582
7fd59977 1583 /* Local variables */
1ef32e96
RL
1584 doublereal alinu, blinu, alinv, blinv;
1585 integer ii, ibb;
7fd59977 1586
1587/* ***********************************************************************
1588 */
1589
0d969553 1590/* FUNCTION : */
7fd59977 1591/* ---------- */
0d969553
Y
1592/* Normalization of parameters of an iso, starting from */
1593/* parametric block and parameters on (-1,1). */
7fd59977 1594
0d969553 1595/* KEYWORDS : */
7fd59977 1596/* ----------- */
1597/* TOUS,AB_SPECIFI :: ISO&,POINT&,NORMALISATION,&POINT,&ISO */
1598
0d969553 1599/* INPUT ARGUMENTS : */
7fd59977 1600/* -------------------- */
0d969553
Y
1601/* NBROOT: Nb of points of discretisation INSIDE the iso */
1602/* defined on (-1,1). */
1603/* ROOTLG: Table of discretization parameters on )-1,1( */
1604/* of the iso. */
1605/* UVFONC: Block of definition of the iso */
1606/* ISOFAV: = 1, this is iso-u; =2, this is iso-v. */
1607
1608/* OUTPUT ARGUMENTS : */
7fd59977 1609/* --------------------- */
0d969553 1610/* TTABLE: Table of parameters renormalized on UVFONC of the iso.
7fd59977 1611*/
1612/* IERCOD: = 0, OK */
0d969553 1613/* = 1, ISOFAV is out of allowed values. */
7fd59977 1614
7fd59977 1615/* > */
1616/* **********************************************************************
1617*/
0d969553 1618/* Name of the routine */
7fd59977 1619
1620
1621 /* Parameter adjustments */
1622 --rootlg;
1623 uvfonc -= 3;
1624
1625 /* Function Body */
1626 ibb = AdvApp2Var_SysBase::mnfndeb_();
1627 if (ibb >= 3) {
1628 AdvApp2Var_SysBase::mgenmsg_("MMA1NOP", 7L);
1629 }
1630
1631 alinu = (uvfonc[4] - uvfonc[3]) / 2.;
1632 blinu = (uvfonc[4] + uvfonc[3]) / 2.;
1633 alinv = (uvfonc[6] - uvfonc[5]) / 2.;
1634 blinv = (uvfonc[6] + uvfonc[5]) / 2.;
1635
1636 if (*isofav == 1) {
1637 ttable[0] = uvfonc[5];
1638 i__1 = *nbroot;
1639 for (ii = 1; ii <= i__1; ++ii) {
1640 ttable[ii] = alinv * rootlg[ii] + blinv;
1641/* L100: */
1642 }
1643 ttable[*nbroot + 1] = uvfonc[6];
1644 } else if (*isofav == 2) {
1645 ttable[0] = uvfonc[3];
1646 i__1 = *nbroot;
1647 for (ii = 1; ii <= i__1; ++ii) {
1648 ttable[ii] = alinu * rootlg[ii] + blinu;
1649/* L200: */
1650 }
1651 ttable[*nbroot + 1] = uvfonc[4];
1652 } else {
1653 goto L9100;
1654 }
1655
1656 goto L9999;
1657
1658/* ------------------------------ THE END -------------------------------
1659*/
1660
1661L9100:
1662 *iercod = 1;
1663 goto L9999;
1664
1665L9999:
1666 if (*iercod != 0) {
1667 AdvApp2Var_SysBase::maermsg_("MMA1NOP", iercod, 7L);
1668 }
1669 if (ibb >= 3) {
1670 AdvApp2Var_SysBase::mgsomsg_("MMA1NOP", 7L);
1671 }
1672
1673 return 0 ;
1674
1675} /* mma1nop_ */
1676
1677//=======================================================================
1678//function : mma2ac1_
1679//purpose :
1680//=======================================================================
1681int AdvApp2Var_ApproxF2var::mma2ac1_(integer const *ndimen,
1682 integer const *mxujac,
1683 integer const *mxvjac,
1684 integer const *iordru,
1685 integer const *iordrv,
1686 doublereal const *contr1,
1687 doublereal const * contr2,
1688 doublereal const *contr3,
1689 doublereal const *contr4,
1690 doublereal const *uhermt,
1691 doublereal const *vhermt,
1692 doublereal *patjac)
1693
1694{
1695 /* System generated locals */
1696 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
1697 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
1698 contr4_dim1, contr4_dim2, contr4_offset, uhermt_dim1,
1699 uhermt_offset, vhermt_dim1, vhermt_offset, patjac_dim1,
1700 patjac_dim2, patjac_offset, i__1, i__2, i__3, i__4, i__5;
41194117 1701
7fd59977 1702 /* Local variables */
1ef32e96
RL
1703 logical ldbg;
1704 integer ndgu, ndgv;
1705 doublereal bidu1, bidu2, bidv1, bidv2;
1706 integer ioru1, iorv1, ii, nd, jj, ku, kv;
1707 doublereal cnt1, cnt2, cnt3, cnt4;
7fd59977 1708
1709/* **********************************************************************
1710*/
1711
0d969553 1712/* FUNCTION : */
7fd59977 1713/* ---------- */
0d969553 1714/* Add polynoms of edge constraints. */
7fd59977 1715
0d969553 1716/* KEYWORDS : */
7fd59977 1717/* ----------- */
1718/* TOUS,AB_SPECIFI::POINT&,CONTRAINTE&,ADDITION,&POLYNOME */
1719
0d969553 1720/* INPUT ARGUMENTS : */
7fd59977 1721/* ------------------ */
0d969553
Y
1722/* NDIMEN: Dimension of the space. */
1723/* MXUJAC: Max degree of the polynom of approximation by U. The */
1724/* representation in the orthogonal base starts from degree */
1725/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1726/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1727/* MXVJAC: Max degree of the polynom of approximation by V. The */
1728/* representation in the orthogonal base starts from degree */
1729/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1730/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1731/* IORDRU: Order of the base of Jacobi (-1,0,1 or 2) by U. Corresponds */
1732/* to the step of constraints: C0, C1 or C2. */
1733/* IORDRV: Order of the base of Jacobi (-1,0,1 or 2) by V. Corresponds */
1734/* to the step of constraints: C0, C1 or C2. */
1735/* CONTR1: Contains, if IORDRU and IORDRV>=0, the values at the */
1736/* extremities of F(U0,V0) and its derivatives. */
1737/* CONTR2: Contains, if IORDRU and IORDRV>=0, the values at the */
1738/* extremities of F(U1,V0) and its derivatives. */
1739/* CONTR3: Contains, if IORDRU and IORDRV>=0, the values at the */
1740/* extremities of F(U0,V1) and its derivatives. */
1741/* CONTR4: Contains, if IORDRU and IORDRV>=0, the values at the */
1742/* extremities of F(U1,V1) and its derivatives. */
1743/* UHERMT: Coeff. of Hermit polynoms of order IORDRU. */
1744/* VHERMT: Coeff. of Hermit polynoms of order IORDRV. */
1745/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
1746/* of F(u,v) WITHOUT taking into account the constraints. */
1747
1748/* OUTPUT ARGUMENTS : */
7fd59977 1749/* ------------------- */
0d969553
Y
1750/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
1751/* of F(u,v) WITH taking into account of constraints. */
7fd59977 1752/* > */
1753/* **********************************************************************
1754*/
0d969553 1755/* Name of the routine */
7fd59977 1756
0d969553 1757/* --------------------------- Initialization --------------------------
7fd59977 1758*/
1759
1760 /* Parameter adjustments */
1761 patjac_dim1 = *mxujac + 1;
1762 patjac_dim2 = *mxvjac + 1;
1763 patjac_offset = patjac_dim1 * patjac_dim2;
1764 patjac -= patjac_offset;
1765 uhermt_dim1 = (*iordru << 1) + 2;
1766 uhermt_offset = uhermt_dim1;
1767 uhermt -= uhermt_offset;
1768 vhermt_dim1 = (*iordrv << 1) + 2;
1769 vhermt_offset = vhermt_dim1;
1770 vhermt -= vhermt_offset;
1771 contr4_dim1 = *ndimen;
1772 contr4_dim2 = *iordru + 2;
1773 contr4_offset = contr4_dim1 * (contr4_dim2 + 1) + 1;
1774 contr4 -= contr4_offset;
1775 contr3_dim1 = *ndimen;
1776 contr3_dim2 = *iordru + 2;
1777 contr3_offset = contr3_dim1 * (contr3_dim2 + 1) + 1;
1778 contr3 -= contr3_offset;
1779 contr2_dim1 = *ndimen;
1780 contr2_dim2 = *iordru + 2;
1781 contr2_offset = contr2_dim1 * (contr2_dim2 + 1) + 1;
1782 contr2 -= contr2_offset;
1783 contr1_dim1 = *ndimen;
1784 contr1_dim2 = *iordru + 2;
1785 contr1_offset = contr1_dim1 * (contr1_dim2 + 1) + 1;
1786 contr1 -= contr1_offset;
1787
1788 /* Function Body */
1789 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
1790 if (ldbg) {
1791 AdvApp2Var_SysBase::mgenmsg_("MMA2AC1", 7L);
1792 }
1793
0d969553 1794/* ------------ SUBTRACTION OF ANGULAR CONSTRAINTS -------------------
7fd59977 1795*/
1796
1797 ioru1 = *iordru + 1;
1798 iorv1 = *iordrv + 1;
1799 ndgu = (*iordru << 1) + 1;
1800 ndgv = (*iordrv << 1) + 1;
1801
1802 i__1 = iorv1;
1803 for (jj = 1; jj <= i__1; ++jj) {
1804 i__2 = ioru1;
1805 for (ii = 1; ii <= i__2; ++ii) {
1806 i__3 = *ndimen;
1807 for (nd = 1; nd <= i__3; ++nd) {
1808 cnt1 = contr1[nd + (ii + jj * contr1_dim2) * contr1_dim1];
1809 cnt2 = contr2[nd + (ii + jj * contr2_dim2) * contr2_dim1];
1810 cnt3 = contr3[nd + (ii + jj * contr3_dim2) * contr3_dim1];
1811 cnt4 = contr4[nd + (ii + jj * contr4_dim2) * contr4_dim1];
1812 i__4 = ndgv;
1813 for (kv = 0; kv <= i__4; ++kv) {
1814 bidv1 = vhermt[kv + ((jj << 1) - 1) * vhermt_dim1];
1815 bidv2 = vhermt[kv + (jj << 1) * vhermt_dim1];
1816 i__5 = ndgu;
1817 for (ku = 0; ku <= i__5; ++ku) {
1818 bidu1 = uhermt[ku + ((ii << 1) - 1) * uhermt_dim1];
1819 bidu2 = uhermt[ku + (ii << 1) * uhermt_dim1];
1820 patjac[ku + (kv + nd * patjac_dim2) * patjac_dim1] =
1821 patjac[ku + (kv + nd * patjac_dim2) *
1822 patjac_dim1] - bidu1 * bidv1 * cnt1 - bidu2 *
1823 bidv1 * cnt2 - bidu1 * bidv2 * cnt3 - bidu2 *
1824 bidv2 * cnt4;
1825/* L500: */
1826 }
1827/* L400: */
1828 }
1829/* L300: */
1830 }
1831/* L200: */
1832 }
1833/* L100: */
1834 }
1835
1836/* ------------------------------ The end -------------------------------
1837*/
1838
1839 if (ldbg) {
1840 AdvApp2Var_SysBase::mgsomsg_("MMA2AC1", 7L);
1841 }
1842 return 0;
1843} /* mma2ac1_ */
1844
1845//=======================================================================
1846//function : mma2ac2_
1847//purpose :
1848//=======================================================================
1849int AdvApp2Var_ApproxF2var::mma2ac2_(const integer *ndimen,
1850 const integer *mxujac,
1851 const integer *mxvjac,
1852 const integer *iordrv,
1853 const integer *nclimu,
1854 const integer *ncfiv1,
1855 const doublereal *crbiv1,
1856 const integer *ncfiv2,
1857 const doublereal *crbiv2,
1858 const doublereal *vhermt,
1859 doublereal *patjac)
1860
1861{
1862 /* System generated locals */
1863 integer crbiv1_dim1, crbiv1_dim2, crbiv1_offset, crbiv2_dim1, crbiv2_dim2,
1864 crbiv2_offset, patjac_dim1, patjac_dim2, patjac_offset,
1865 vhermt_dim1, vhermt_offset, i__1, i__2, i__3, i__4;
41194117 1866
7fd59977 1867 /* Local variables */
1ef32e96
RL
1868 logical ldbg;
1869 integer ndgv1, ndgv2, ii, jj, nd, kk;
1870 doublereal bid1, bid2;
7fd59977 1871
1872/* **********************************************************************
1873*/
1874
0d969553 1875/* FUNCTION : */
7fd59977 1876/* ---------- */
0d969553 1877/* Add polynoms of constraints */
7fd59977 1878
0d969553 1879/* KEYWORDS : */
7fd59977 1880/* ----------- */
0d969553 1881/* FUNCTION,APPROXIMATION,COEFFICIENT,POLYNOM */
7fd59977 1882
0d969553 1883/* INPUT ARGUMENTS : */
7fd59977 1884/* ------------------ */
0d969553
Y
1885/* NDIMEN: Dimension of the space. */
1886/* MXUJAC: Max degree of the polynom of approximation by U. The */
1887/* representation in the orthogonal base starts from degree */
1888/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1889/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1890/* MXVJAC: Max degree of the polynom of approximation by V. The */
1891/* representation in the orthogonal base starts from degree */
1892/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1893/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1894/* IORDRV: Order of the base of Jacobi (-1,0,1 or 2) by V. Corresponds */
1895/* to the step of constraints: C0, C1 or C2. */
1896/* NCLIMU LIMIT nb of coeff by u of the solution P(u,v)
1897* NCFIV1: Nb of Coeff. of curves stored in CRBIV1. */
1898/* CRBIV1: Table of coeffs of the approximation of iso-V0 and its */
1899/* derivatives till order IORDRV. */
1900/* NCFIV2: Nb of Coeff. of curves stored in CRBIV2. */
1901/* CRBIV2: Table of coeffs of approximation of iso-V1 and its */
1902/* derivatives till order IORDRV. */
1903/* VHERMT: Coeff. of Hermit polynoms of order IORDRV. */
1904/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
1905/* of F(u,v) WITHOUT taking into account the constraints. */
1906
1907/* OUTPUT ARGUMENTS : */
7fd59977 1908/* ------------------- */
0d969553
Y
1909/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
1910/* of F(u,v) WITH taking into account of constraints. */
1911/* > *//*
7fd59977 1912
7fd59977 1913
7fd59977 1914/* > */
1915/* **********************************************************************
1916*/
0d969553 1917/* Name of the routine */
7fd59977 1918
1919/* --------------------------- Initialisations --------------------------
1920*/
1921
1922 /* Parameter adjustments */
1923 patjac_dim1 = *mxujac + 1;
1924 patjac_dim2 = *mxvjac + 1;
1925 patjac_offset = patjac_dim1 * patjac_dim2;
1926 patjac -= patjac_offset;
1927 vhermt_dim1 = (*iordrv << 1) + 2;
1928 vhermt_offset = vhermt_dim1;
1929 vhermt -= vhermt_offset;
1930 --ncfiv2;
1931 --ncfiv1;
1932 crbiv2_dim1 = *nclimu;
1933 crbiv2_dim2 = *ndimen;
1934 crbiv2_offset = crbiv2_dim1 * (crbiv2_dim2 + 1);
1935 crbiv2 -= crbiv2_offset;
1936 crbiv1_dim1 = *nclimu;
1937 crbiv1_dim2 = *ndimen;
1938 crbiv1_offset = crbiv1_dim1 * (crbiv1_dim2 + 1);
1939 crbiv1 -= crbiv1_offset;
1940
1941 /* Function Body */
1942 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
1943 if (ldbg) {
1944 AdvApp2Var_SysBase::mgenmsg_("MMA2AC2", 7L);
1945 }
1946
0d969553 1947/* ------------ ADDING of coeff by u of curves, by v of Hermit --------
7fd59977 1948*/
1949
1950 i__1 = *iordrv + 1;
1951 for (ii = 1; ii <= i__1; ++ii) {
1952 ndgv1 = ncfiv1[ii] - 1;
1953 ndgv2 = ncfiv2[ii] - 1;
1954 i__2 = *ndimen;
1955 for (nd = 1; nd <= i__2; ++nd) {
1956 i__3 = (*iordrv << 1) + 1;
1957 for (jj = 0; jj <= i__3; ++jj) {
1958 bid1 = vhermt[jj + ((ii << 1) - 1) * vhermt_dim1];
1959 i__4 = ndgv1;
1960 for (kk = 0; kk <= i__4; ++kk) {
1961 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
1962 bid1 * crbiv1[kk + (nd + ii * crbiv1_dim2) *
1963 crbiv1_dim1];
1964/* L400: */
1965 }
1966 bid2 = vhermt[jj + (ii << 1) * vhermt_dim1];
1967 i__4 = ndgv2;
1968 for (kk = 0; kk <= i__4; ++kk) {
1969 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
1970 bid2 * crbiv2[kk + (nd + ii * crbiv2_dim2) *
1971 crbiv2_dim1];
1972/* L500: */
1973 }
1974/* L300: */
1975 }
1976/* L200: */
1977 }
1978/* L100: */
1979 }
1980
1981/* ------------------------------ The end -------------------------------
1982*/
1983
1984 if (ldbg) {
1985 AdvApp2Var_SysBase::mgsomsg_("MMA2AC2", 7L);
1986 }
1987 return 0;
1988} /* mma2ac2_ */
1989
1990
1991//=======================================================================
1992//function : mma2ac3_
1993//purpose :
1994//=======================================================================
1995int AdvApp2Var_ApproxF2var::mma2ac3_(const integer *ndimen,
1996 const integer *mxujac,
1997 const integer *mxvjac,
1998 const integer *iordru,
1999 const integer *nclimv,
2000 const integer *ncfiu1,
2001 const doublereal * crbiu1,
2002 const integer *ncfiu2,
2003 const doublereal *crbiu2,
2004 const doublereal *uhermt,
2005 doublereal *patjac)
2006
2007{
2008 /* System generated locals */
2009 integer crbiu1_dim1, crbiu1_dim2, crbiu1_offset, crbiu2_dim1, crbiu2_dim2,
2010 crbiu2_offset, patjac_dim1, patjac_dim2, patjac_offset,
2011 uhermt_dim1, uhermt_offset, i__1, i__2, i__3, i__4;
41194117 2012
7fd59977 2013 /* Local variables */
1ef32e96
RL
2014 logical ldbg;
2015 integer ndgu1, ndgu2, ii, jj, nd, kk;
2016 doublereal bid1, bid2;
7fd59977 2017
2018/* **********************************************************************
2019*/
2020
0d969553 2021/* FUNCTION : */
7fd59977 2022/* ---------- */
2023/* Ajout des polynomes de contraintes */
2024
0d969553 2025/* KEYWORDS : */
7fd59977 2026/* ----------- */
2027/* FONCTION,APPROXIMATION,COEFFICIENT,POLYNOME */
2028
0d969553 2029/* INPUT ARGUMENTS : */
7fd59977 2030/* ------------------ */
0d969553
Y
2031/* NDIMEN: Dimension of the space. */
2032/* MXUJAC: Max degree of the polynom of approximation by U. The */
2033/* representation in the orthogonal base starts from degree */
2034/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
2035/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
2036/* MXVJAC: Max degree of the polynom of approximation by V. The */
2037/* representation in the orthogonal base starts from degree */
2038/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
2039/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
2040/* IORDRU: Order of the base of Jacobi (-1,0,1 or 2) by U. Corresponds */
2041/* to the step of constraints: C0, C1 or C2. */
2042/* NCLIMV LIMIT nb of coeff by v of the solution P(u,v)
2043* NCFIU1: Nb of Coeff. of curves stored in CRBIU1. */
2044/* CRBIU1: Table of coeffs of the approximation of iso-U0 and its */
2045/* derivatives till order IORDRU. */
2046/* NCFIU2: Nb of Coeff. of curves stored in CRBIU2. */
2047/* CRBIU2: Table of coeffs of approximation of iso-U1 and its */
2048/* derivatives till order IORDRU */
2049/* UHERMT: Coeff. of Hermit polynoms of order IORDRU. */
2050/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
2051/* of F(u,v) WITHOUT taking into account the constraints. */
2052
2053/* OUTPUT ARGUMENTS : */
7fd59977 2054/* ------------------- */
0d969553
Y
2055/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
2056/* of F(u,v) WITH taking into account of constraints. */
7fd59977 2057
7fd59977 2058
7fd59977 2059/* > */
2060/* **********************************************************************
2061*/
0d969553 2062/* The name of the routine */
7fd59977 2063
0d969553 2064/* --------------------------- Initializations --------------------------
7fd59977 2065*/
2066
2067 /* Parameter adjustments */
2068 patjac_dim1 = *mxujac + 1;
2069 patjac_dim2 = *mxvjac + 1;
2070 patjac_offset = patjac_dim1 * patjac_dim2;
2071 patjac -= patjac_offset;
2072 uhermt_dim1 = (*iordru << 1) + 2;
2073 uhermt_offset = uhermt_dim1;
2074 uhermt -= uhermt_offset;
2075 --ncfiu2;
2076 --ncfiu1;
2077 crbiu2_dim1 = *nclimv;
2078 crbiu2_dim2 = *ndimen;
2079 crbiu2_offset = crbiu2_dim1 * (crbiu2_dim2 + 1);
2080 crbiu2 -= crbiu2_offset;
2081 crbiu1_dim1 = *nclimv;
2082 crbiu1_dim2 = *ndimen;
2083 crbiu1_offset = crbiu1_dim1 * (crbiu1_dim2 + 1);
2084 crbiu1 -= crbiu1_offset;
2085
2086 /* Function Body */
2087 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
2088 if (ldbg) {
2089 AdvApp2Var_SysBase::mgenmsg_("MMA2AC3", 7L);
2090 }
2091
0d969553 2092/* ------------ ADDING of coeff by u of curves, by v of Hermit --------
7fd59977 2093*/
2094
2095 i__1 = *iordru + 1;
2096 for (ii = 1; ii <= i__1; ++ii) {
2097 ndgu1 = ncfiu1[ii] - 1;
2098 ndgu2 = ncfiu2[ii] - 1;
2099 i__2 = *ndimen;
2100 for (nd = 1; nd <= i__2; ++nd) {
2101 i__3 = ndgu1;
2102 for (jj = 0; jj <= i__3; ++jj) {
2103 bid1 = crbiu1[jj + (nd + ii * crbiu1_dim2) * crbiu1_dim1];
2104 i__4 = (*iordru << 1) + 1;
2105 for (kk = 0; kk <= i__4; ++kk) {
2106 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
2107 bid1 * uhermt[kk + ((ii << 1) - 1) * uhermt_dim1];
2108/* L400: */
2109 }
2110/* L300: */
2111 }
2112 i__3 = ndgu2;
2113 for (jj = 0; jj <= i__3; ++jj) {
2114 bid2 = crbiu2[jj + (nd + ii * crbiu2_dim2) * crbiu2_dim1];
2115 i__4 = (*iordru << 1) + 1;
2116 for (kk = 0; kk <= i__4; ++kk) {
2117 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
2118 bid2 * uhermt[kk + (ii << 1) * uhermt_dim1];
2119/* L600: */
2120 }
2121/* L500: */
2122 }
2123
2124/* L200: */
2125 }
2126/* L100: */
2127 }
2128
2129/* ------------------------------ The end -------------------------------
2130*/
2131
2132 if (ldbg) {
2133 AdvApp2Var_SysBase::mgsomsg_("MMA2AC3", 7L);
2134 }
2135 return 0;
2136} /* mma2ac3_ */
2137
2138//=======================================================================
2139//function : mma2can_
2140//purpose :
2141//=======================================================================
2142int AdvApp2Var_ApproxF2var::mma2can_(const integer *ncfmxu,
2143 const integer *ncfmxv,
2144 const integer *ndimen,
2145 const integer *iordru,
2146 const integer *iordrv,
2147 const integer *ncoefu,
2148 const integer *ncoefv,
2149 const doublereal *patjac,
2150 doublereal *pataux,
2151 doublereal *patcan,
2152 integer *iercod)
2153
2154{
2155 /* System generated locals */
2156 integer patjac_dim1, patjac_dim2, patjac_offset, patcan_dim1, patcan_dim2,
2157 patcan_offset, i__1, i__2;
41194117 2158
7fd59977 2159 /* Local variables */
1ef32e96
RL
2160 logical ldbg;
2161 integer ilon1, ilon2, ii, nd;
7fd59977 2162
2163/* **********************************************************************
2164*/
2165
0d969553 2166/* FUNCTION : */
7fd59977 2167/* ---------- */
0d969553
Y
2168/* Change of Jacobi base to canonical (-1,1) and writing in a greater */
2169/* table. */
7fd59977 2170
0d969553 2171/* KEYWORDS : */
7fd59977 2172/* ----------- */
0d969553 2173/* ALL,AB_SPECIFI,CARREAU&,CONVERSION,JACOBI,CANNONIQUE,&CARREAU */
7fd59977 2174
0d969553 2175/* INPUT ARGUMENTS : */
7fd59977 2176/* -------------------- */
0d969553
Y
2177/* NCFMXU: Dimension by U of resulting table PATCAN */
2178/* NCFMXV: Dimension by V of resulting table PATCAN */
2179/* NDIMEN: Dimension of the workspace. */
2180/* IORDRU: Order of constraint by U */
2181/* IORDRV: Order of constraint by V. */
2182/* NCOEFU: Nb of coeff by U of square PATJAC */
2183/* NCOEFV: Nb of coeff by V of square PATJAC */
2184/* PATJAC: Square in the base of Jacobi of order IORDRU by U and */
2185/* IORDRV by V. */
2186
2187/* OUTPUT ARGUMENTS : */
7fd59977 2188/* --------------------- */
0d969553
Y
2189/* PATAUX: Auxiliary Table. */
2190/* PATCAN: Table of coefficients in the canonic base. */
2191/* IERCOD: Error code. */
2192/* = 0, everything goes well, and all things are equal. */
2193/* = 1, the program refuses to process with incorrect input arguments */
2194
2195
2196/* COMMONS USED : */
7fd59977 2197/* ------------------ */
2198
0d969553 2199/* REFERENCES CALLED : */
7fd59977 2200/* --------------------- */
2201
0d969553 2202/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2203/* ----------------------------------- */
7fd59977 2204/* > */
2205/* **********************************************************************
2206*/
2207
2208
2209 /* Parameter adjustments */
2210 patcan_dim1 = *ncfmxu;
2211 patcan_dim2 = *ncfmxv;
2212 patcan_offset = patcan_dim1 * (patcan_dim2 + 1) + 1;
2213 patcan -= patcan_offset;
2214 --pataux;
2215 patjac_dim1 = *ncoefu;
2216 patjac_dim2 = *ncoefv;
2217 patjac_offset = patjac_dim1 * (patjac_dim2 + 1) + 1;
2218 patjac -= patjac_offset;
2219
2220 /* Function Body */
2221 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 2;
2222 if (ldbg) {
2223 AdvApp2Var_SysBase::mgenmsg_("MMA2CAN", 7L);
2224 }
2225 *iercod = 0;
2226
2227 if (*iordru < -1 || *iordru > 2) {
2228 goto L9100;
2229 }
2230 if (*iordrv < -1 || *iordrv > 2) {
2231 goto L9100;
2232 }
2233 if (*ncoefu > *ncfmxu || *ncoefv > *ncfmxv) {
2234 goto L9100;
2235 }
2236
0d969553 2237/* --> Pass to canonic base (-1,1) */
7fd59977 2238 mmjacpt_(ndimen, ncoefu, ncoefv, iordru, iordrv, &patjac[patjac_offset], &
2239 pataux[1], &patcan[patcan_offset]);
2240
0d969553 2241/* --> Write all in a greater table */
fadcea2c
RL
2242 AdvApp2Var_MathBase::mmfmca8_(ncoefu,
2243 ncoefv,
2244 ndimen,
2245 ncfmxu,
2246 ncfmxv,
2247 ndimen,
2248 &patcan[patcan_offset],
2249 &patcan[patcan_offset]);
7fd59977 2250
0d969553 2251/* --> Complete with zeros the resulting table. */
7fd59977 2252 ilon1 = *ncfmxu - *ncoefu;
2253 ilon2 = *ncfmxu * (*ncfmxv - *ncoefv);
2254 i__1 = *ndimen;
2255 for (nd = 1; nd <= i__1; ++nd) {
2256 if (ilon1 > 0) {
2257 i__2 = *ncoefv;
2258 for (ii = 1; ii <= i__2; ++ii) {
2259 AdvApp2Var_SysBase::mvriraz_(&ilon1,
fadcea2c 2260 &patcan[*ncoefu + 1 + (ii + nd * patcan_dim2) * patcan_dim1]);
7fd59977 2261/* L110: */
2262 }
2263 }
2264 if (ilon2 > 0) {
2265 AdvApp2Var_SysBase::mvriraz_(&ilon2,
fadcea2c 2266 &patcan[(*ncoefv + 1 + nd * patcan_dim2) * patcan_dim1 + 1]);
7fd59977 2267 }
2268/* L100: */
2269 }
2270
2271 goto L9999;
2272
0d969553 2273/* ----------------------
7fd59977 2274*/
2275
2276L9100:
2277 *iercod = 1;
2278 goto L9999;
2279
2280L9999:
2281 AdvApp2Var_SysBase::maermsg_("MMA2CAN", iercod, 7L);
2282 if (ldbg) {
2283 AdvApp2Var_SysBase::mgsomsg_("MMA2CAN", 7L);
2284 }
2285 return 0 ;
2286} /* mma2can_ */
2287
2288//=======================================================================
2289//function : mma2cd1_
2290//purpose :
2291//=======================================================================
2292int mma2cd1_(integer *ndimen,
2293 integer *nbpntu,
2294 doublereal *urootl,
2295 integer *nbpntv,
2296 doublereal *vrootl,
2297 integer *iordru,
2298 integer *iordrv,
2299 doublereal *contr1,
2300 doublereal *contr2,
2301 doublereal *contr3,
2302 doublereal *contr4,
2303 doublereal *fpntbu,
2304 doublereal *fpntbv,
2305 doublereal *uhermt,
2306 doublereal *vhermt,
2307 doublereal *sosotb,
2308 doublereal *soditb,
2309 doublereal *disotb,
2310 doublereal *diditb)
2311
2312{
1ef32e96 2313 integer c__1 = 1;
41194117 2314
7fd59977 2315/* System generated locals */
2316 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
2317 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
2318 contr4_dim1, contr4_dim2, contr4_offset, uhermt_dim1,
2319 uhermt_offset, vhermt_dim1, vhermt_offset, fpntbu_dim1,
2320 fpntbu_offset, fpntbv_dim1, fpntbv_offset, sosotb_dim1,
2321 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2322 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2323 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4,
2324 i__5;
2325
2326 /* Local variables */
1ef32e96 2327 integer ncfhu, ncfhv, nuroo, nvroo, nd, ii, jj, kk, ll, ibb, kkm,
7fd59977 2328 llm, kkp, llp;
1ef32e96
RL
2329 doublereal bid1, bid2, bid3, bid4;
2330 doublereal diu1, diu2, div1, div2, sou1, sou2, sov1, sov2;
7fd59977 2331
7fd59977 2332/* **********************************************************************
2333*/
2334
0d969553 2335/* FUNCTION : */
7fd59977 2336/* ---------- */
0d969553
Y
2337/* Discretisation on the parameters of polynoms of interpolation */
2338/* of constraints at the corners of order IORDRE. */
7fd59977 2339
0d969553 2340/* KEYWORDS : */
7fd59977 2341/* ----------- */
2342/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
2343
0d969553 2344/* INPUT ARGUMENTS : */
7fd59977 2345/* ------------------ */
0d969553
Y
2346/* NDIMEN: Dimension of the space. */
2347/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
2348/* This is also the nb of root of Legendre polynom where discretization is done. */
2349/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
2350*/
2351/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
2352/* This is also the nb of root of Legendre polynom where discretization is done. */
2353/* VROOTL: Table of discretization parameters on (-1,1) by V.
2354/* IORDRU: Order of constraint imposed at the extremities of iso-V */
2355/* = 0, calculate the extremities of iso-V */
2356/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
2357/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
2358/* IORDRV: Order of constraint imposed at the extremities of iso-U */
2359/* = 0, calculate the extremities of iso-U */
2360/* = 1, calculate, additionally, the 1st derivative in the direction of iso-U */
2361/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-U */
2362/* CONTR1: Contains, if IORDRU and IORDRV>=0, the values at the */
2363/* extremities of F(U0,V0) and its derivatives. */
2364/* CONTR2: Contains, if IORDRU and IORDRV>=0, the values at the */
2365/* extremities of F(U1,V0) and its derivatives. */
2366/* CONTR3: Contains, if IORDRU and IORDRV>=0, the values at the */
2367/* extremities of F(U0,V1) and its derivatives. */
2368/* CONTR4: Contains, if IORDRU and IORDRV>=0, the values at the */
2369/* extremities of F(U1,V1) and its derivatives. */
2370/* SOSOTB: Preinitialized table (input/output argument). */
2371/* DISOTB: Preinitialized table (input/output argument). */
2372/* SODITB: Preinitialized table (input/output argument). */
2373/* DIDITB: Preinitialized table (input/output argument) */
2374
2375/* OUTPUT ARGUMENTS : */
7fd59977 2376/* ------------------- */
0d969553
Y
2377/* FPNTBU: Auxiliary table. */
2378/* FPNTBV: Auxiliary table. */
2379/* UHERMT: Table of 2*(IORDRU+1) coeff. of 2*(IORDRU+1) polynoms of Hermite. */
2380/* VHERMT: Table of 2*(IORDRV+1) coeff. of 2*(IORDRV+1) polynoms of Hermite. */
2381/* SOSOTB: Table where the terms of constraints are added */
7fd59977 2382/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2383/* with ui and vj positive roots of the Legendre polynom */
2384/* of degree NBPNTU and NBPNTV respectively. */
2385/* DISOTB: Table where the terms of constraints are added */
7fd59977 2386/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2387/* with ui and vj positive roots of the polynom of Legendre */
2388/* of degree NBPNTU and NBPNTV respectively. */
2389/* SODITB: Table where the terms of constraints are added */
7fd59977 2390/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2391/* with ui and vj positive roots of the polynom of Legendre */
2392/* of degree NBPNTU and NBPNTV respectively. */
2393/* DIDITB: Table where the terms of constraints are added */
7fd59977 2394/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2395/* with ui and vj positive roots of the polynom of Legendre */
2396/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 2397
0d969553 2398/* COMMONS USED : */
7fd59977 2399/* ---------------- */
2400
0d969553 2401/* REFERENCES CALLED : */
7fd59977 2402/* ----------------------- */
2403
0d969553 2404/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2405/* ----------------------------------- */
2406
7fd59977 2407/* > */
2408/* **********************************************************************
2409*/
2410
0d969553 2411/* Name of the routine */
7fd59977 2412
2413
2414 /* Parameter adjustments */
2415 --urootl;
2416 diditb_dim1 = *nbpntu / 2 + 1;
2417 diditb_dim2 = *nbpntv / 2 + 1;
2418 diditb_offset = diditb_dim1 * diditb_dim2;
2419 diditb -= diditb_offset;
2420 disotb_dim1 = *nbpntu / 2;
2421 disotb_dim2 = *nbpntv / 2;
2422 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
2423 disotb -= disotb_offset;
2424 soditb_dim1 = *nbpntu / 2;
2425 soditb_dim2 = *nbpntv / 2;
2426 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
2427 soditb -= soditb_offset;
2428 sosotb_dim1 = *nbpntu / 2 + 1;
2429 sosotb_dim2 = *nbpntv / 2 + 1;
2430 sosotb_offset = sosotb_dim1 * sosotb_dim2;
2431 sosotb -= sosotb_offset;
2432 --vrootl;
2433 uhermt_dim1 = (*iordru << 1) + 2;
2434 uhermt_offset = uhermt_dim1;
2435 uhermt -= uhermt_offset;
2436 fpntbu_dim1 = *nbpntu;
2437 fpntbu_offset = fpntbu_dim1 + 1;
2438 fpntbu -= fpntbu_offset;
2439 vhermt_dim1 = (*iordrv << 1) + 2;
2440 vhermt_offset = vhermt_dim1;
2441 vhermt -= vhermt_offset;
2442 fpntbv_dim1 = *nbpntv;
2443 fpntbv_offset = fpntbv_dim1 + 1;
2444 fpntbv -= fpntbv_offset;
2445 contr4_dim1 = *ndimen;
2446 contr4_dim2 = *iordru + 2;
2447 contr4_offset = contr4_dim1 * (contr4_dim2 + 1) + 1;
2448 contr4 -= contr4_offset;
2449 contr3_dim1 = *ndimen;
2450 contr3_dim2 = *iordru + 2;
2451 contr3_offset = contr3_dim1 * (contr3_dim2 + 1) + 1;
2452 contr3 -= contr3_offset;
2453 contr2_dim1 = *ndimen;
2454 contr2_dim2 = *iordru + 2;
2455 contr2_offset = contr2_dim1 * (contr2_dim2 + 1) + 1;
2456 contr2 -= contr2_offset;
2457 contr1_dim1 = *ndimen;
2458 contr1_dim2 = *iordru + 2;
2459 contr1_offset = contr1_dim1 * (contr1_dim2 + 1) + 1;
2460 contr1 -= contr1_offset;
2461
2462 /* Function Body */
2463 ibb = AdvApp2Var_SysBase::mnfndeb_();
2464 if (ibb >= 3) {
2465 AdvApp2Var_SysBase::mgenmsg_("MMA2CD1", 7L);
2466 }
2467
0d969553 2468/* ------------------- Discretisation of Hermite polynoms -----------
7fd59977 2469*/
2470
2471 ncfhu = (*iordru + 1) << 1;
2472 i__1 = ncfhu;
2473 for (ii = 1; ii <= i__1; ++ii) {
2474 i__2 = *nbpntu;
2475 for (ll = 1; ll <= i__2; ++ll) {
2476 AdvApp2Var_MathBase::mmmpocur_(&ncfhu, &c__1, &ncfhu, &uhermt[ii * uhermt_dim1], &
2477 urootl[ll], &fpntbu[ll + ii * fpntbu_dim1]);
2478/* L20: */
2479 }
2480/* L10: */
2481 }
2482 ncfhv = (*iordrv + 1) << 1;
2483 i__1 = ncfhv;
2484 for (jj = 1; jj <= i__1; ++jj) {
2485 i__2 = *nbpntv;
2486 for (kk = 1; kk <= i__2; ++kk) {
2487 AdvApp2Var_MathBase::mmmpocur_(&ncfhv, &c__1, &ncfhv, &vhermt[jj * vhermt_dim1], &
2488 vrootl[kk], &fpntbv[kk + jj * fpntbv_dim1]);
2489/* L40: */
2490 }
2491/* L30: */
2492 }
2493
0d969553 2494/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 2495*/
2496
2497 nuroo = *nbpntu / 2;
2498 nvroo = *nbpntv / 2;
2499 i__1 = *ndimen;
2500 for (nd = 1; nd <= i__1; ++nd) {
2501
2502 i__2 = *iordrv + 1;
2503 for (jj = 1; jj <= i__2; ++jj) {
2504 i__3 = *iordru + 1;
2505 for (ii = 1; ii <= i__3; ++ii) {
2506 bid1 = contr1[nd + (ii + jj * contr1_dim2) * contr1_dim1];
2507 bid2 = contr2[nd + (ii + jj * contr2_dim2) * contr2_dim1];
2508 bid3 = contr3[nd + (ii + jj * contr3_dim2) * contr3_dim1];
2509 bid4 = contr4[nd + (ii + jj * contr4_dim2) * contr4_dim1];
2510
2511 i__4 = nvroo;
2512 for (kk = 1; kk <= i__4; ++kk) {
2513 kkp = (*nbpntv + 1) / 2 + kk;
2514 kkm = nvroo - kk + 1;
2515 sov1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] +
2516 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2517 div1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] -
2518 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2519 sov2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] + fpntbv[kkm
2520 + (jj << 1) * fpntbv_dim1];
2521 div2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] - fpntbv[kkm
2522 + (jj << 1) * fpntbv_dim1];
2523 i__5 = nuroo;
2524 for (ll = 1; ll <= i__5; ++ll) {
2525 llp = (*nbpntu + 1) / 2 + ll;
2526 llm = nuroo - ll + 1;
2527 sou1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] +
2528 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2529 diu1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] -
2530 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2531 sou2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] + fpntbu[
2532 llm + (ii << 1) * fpntbu_dim1];
2533 diu2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] - fpntbu[
2534 llm + (ii << 1) * fpntbu_dim1];
2535 sosotb[ll + (kk + nd * sosotb_dim2) * sosotb_dim1] =
2536 sosotb[ll + (kk + nd * sosotb_dim2) *
2537 sosotb_dim1] - bid1 * sou1 * sov1 - bid2 *
2538 sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2539 sou2 * sov2;
2540 soditb[ll + (kk + nd * soditb_dim2) * soditb_dim1] =
2541 soditb[ll + (kk + nd * soditb_dim2) *
2542 soditb_dim1] - bid1 * sou1 * div1 - bid2 *
2543 sou2 * div1 - bid3 * sou1 * div2 - bid4 *
2544 sou2 * div2;
2545 disotb[ll + (kk + nd * disotb_dim2) * disotb_dim1] =
2546 disotb[ll + (kk + nd * disotb_dim2) *
2547 disotb_dim1] - bid1 * diu1 * sov1 - bid2 *
2548 diu2 * sov1 - bid3 * diu1 * sov2 - bid4 *
2549 diu2 * sov2;
2550 diditb[ll + (kk + nd * diditb_dim2) * diditb_dim1] =
2551 diditb[ll + (kk + nd * diditb_dim2) *
2552 diditb_dim1] - bid1 * diu1 * div1 - bid2 *
2553 diu2 * div1 - bid3 * diu1 * div2 - bid4 *
2554 diu2 * div2;
2555/* L450: */
2556 }
2557/* L400: */
2558 }
2559
0d969553
Y
2560/* ------------ Case when the discretization is done only on the roots
2561----------- */
2562/* ---------- of Legendre polynom of uneven degree, 0 is root
7fd59977 2563----------- */
7fd59977 2564
2565 if (*nbpntu % 2 == 1) {
2566 sou1 = fpntbu[nuroo + 1 + ((ii << 1) - 1) * fpntbu_dim1];
2567 sou2 = fpntbu[nuroo + 1 + (ii << 1) * fpntbu_dim1];
2568 i__4 = nvroo;
2569 for (kk = 1; kk <= i__4; ++kk) {
2570 kkp = (*nbpntv + 1) / 2 + kk;
2571 kkm = nvroo - kk + 1;
2572 sov1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] +
2573 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2574 div1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] -
2575 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2576 sov2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] + fpntbv[
2577 kkm + (jj << 1) * fpntbv_dim1];
2578 div2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] - fpntbv[
2579 kkm + (jj << 1) * fpntbv_dim1];
2580 sosotb[(kk + nd * sosotb_dim2) * sosotb_dim1] =
2581 sosotb[(kk + nd * sosotb_dim2) * sosotb_dim1]
2582 - bid1 * sou1 * sov1 - bid2 * sou2 * sov1 -
2583 bid3 * sou1 * sov2 - bid4 * sou2 * sov2;
2584 diditb[(kk + nd * diditb_dim2) * diditb_dim1] =
2585 diditb[(kk + nd * diditb_dim2) * diditb_dim1]
2586 - bid1 * sou1 * div1 - bid2 * sou2 * div1 -
2587 bid3 * sou1 * div2 - bid4 * sou2 * div2;
2588/* L500: */
2589 }
2590 }
2591
2592 if (*nbpntv % 2 == 1) {
2593 sov1 = fpntbv[nvroo + 1 + ((jj << 1) - 1) * fpntbv_dim1];
2594 sov2 = fpntbv[nvroo + 1 + (jj << 1) * fpntbv_dim1];
2595 i__4 = nuroo;
2596 for (ll = 1; ll <= i__4; ++ll) {
2597 llp = (*nbpntu + 1) / 2 + ll;
2598 llm = nuroo - ll + 1;
2599 sou1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] +
2600 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2601 diu1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] -
2602 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2603 sou2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] + fpntbu[
2604 llm + (ii << 1) * fpntbu_dim1];
2605 diu2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] - fpntbu[
2606 llm + (ii << 1) * fpntbu_dim1];
2607 sosotb[ll + nd * sosotb_dim2 * sosotb_dim1] = sosotb[
2608 ll + nd * sosotb_dim2 * sosotb_dim1] - bid1 *
2609 sou1 * sov1 - bid2 * sou2 * sov1 - bid3 *
2610 sou1 * sov2 - bid4 * sou2 * sov2;
2611 diditb[ll + nd * diditb_dim2 * diditb_dim1] = diditb[
2612 ll + nd * diditb_dim2 * diditb_dim1] - bid1 *
2613 diu1 * sov1 - bid2 * diu2 * sov1 - bid3 *
2614 diu1 * sov2 - bid4 * diu2 * sov2;
2615/* L600: */
2616 }
2617 }
2618
2619 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
2620 sou1 = fpntbu[nuroo + 1 + ((ii << 1) - 1) * fpntbu_dim1];
2621 sou2 = fpntbu[nuroo + 1 + (ii << 1) * fpntbu_dim1];
2622 sov1 = fpntbv[nvroo + 1 + ((jj << 1) - 1) * fpntbv_dim1];
2623 sov2 = fpntbv[nvroo + 1 + (jj << 1) * fpntbv_dim1];
2624 sosotb[nd * sosotb_dim2 * sosotb_dim1] = sosotb[nd *
2625 sosotb_dim2 * sosotb_dim1] - bid1 * sou1 * sov1 -
2626 bid2 * sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2627 sou2 * sov2;
2628 diditb[nd * diditb_dim2 * diditb_dim1] = diditb[nd *
2629 diditb_dim2 * diditb_dim1] - bid1 * sou1 * sov1 -
2630 bid2 * sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2631 sou2 * sov2;
2632 }
2633
2634/* L300: */
2635 }
2636/* L200: */
2637 }
2638/* L100: */
2639 }
2640 goto L9999;
2641
2642/* ------------------------------ The End -------------------------------
2643*/
2644
2645L9999:
2646 if (ibb >= 3) {
2647 AdvApp2Var_SysBase::mgsomsg_("MMA2CD1", 7L);
2648 }
2649 return 0;
2650} /* mma2cd1_ */
2651
2652//=======================================================================
2653//function : mma2cd2_
2654//purpose :
2655//=======================================================================
2656int mma2cd2_(integer *ndimen,
2657 integer *nbpntu,
2658 integer *nbpntv,
2659 doublereal *vrootl,
2660 integer *iordrv,
2661 doublereal *sotbv1,
2662 doublereal *sotbv2,
2663 doublereal *ditbv1,
2664 doublereal *ditbv2,
2665 doublereal *fpntab,
2666 doublereal *vhermt,
2667 doublereal *sosotb,
2668 doublereal *soditb,
2669 doublereal *disotb,
2670 doublereal *diditb)
2671
2672{
1ef32e96 2673 integer c__1 = 1;
7fd59977 2674 /* System generated locals */
2675 integer sotbv1_dim1, sotbv1_dim2, sotbv1_offset, sotbv2_dim1, sotbv2_dim2,
2676 sotbv2_offset, ditbv1_dim1, ditbv1_dim2, ditbv1_offset,
2677 ditbv2_dim1, ditbv2_dim2, ditbv2_offset, fpntab_dim1,
2678 fpntab_offset, vhermt_dim1, vhermt_offset, sosotb_dim1,
2679 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2680 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2681 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4;
41194117 2682
7fd59977 2683 /* Local variables */
1ef32e96
RL
2684 integer ncfhv, nuroo, nvroo, ii, nd, jj, kk, ibb, jjm, jjp;
2685 doublereal bid1, bid2, bid3, bid4;
7fd59977 2686
2687/* **********************************************************************
2688*/
0d969553 2689/* FUNCTION : */
7fd59977 2690/* ---------- */
0d969553
Y
2691/* Discretisation on the parameters of polynoms of interpolation */
2692/* of constraints on 2 borders iso-V of order IORDRV. */
7fd59977 2693
0d969553
Y
2694
2695/* KEYWORDS : */
7fd59977 2696/* ----------- */
2697/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
2698
7fd59977 2699
0d969553
Y
2700
2701/* INPUT ARGUMENTS : */
2702/* ------------------ */
2703/* NDIMEN: Dimension of the space. */
2704/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
2705/* This is also the nb of root of Legendre polynom where discretization is done. */
2706/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
2707*/
2708/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
2709/* This is also the nb of root of Legendre polynom where discretization is done. */
2710/* VROOTL: Table of discretization parameters on (-1,1) by V.
2711/* IORDRV: Order of constraint imposed at the extremities of iso-V */
2712/* = 0, calculate the extremities of iso-V */
2713/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
2714/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
2715/* SOTBV1: Table of NBPNTV/2 sums of 2 index points */
2716/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V0. */
2717/* SOTBV2: Table of NBPNTV/2 sums of 2 index points */
2718/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V1. */
2719/* DITBV1: Table of NBPNTV/2 differences of 2 index points */
2720/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V0. */
2721/* DITBV2: Table of NBPNTV/2 differences of 2 index points */
2722/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V1. */
2723/* SOSOTB: Preinitialized table (input/output argument). */
2724/* DISOTB: Preinitialized table (input/output argument). */
2725/* SODITB: Preinitialized table (input/output argument). */
2726/* DIDITB: Preinitialized table (input/output argument) */
2727
2728/* OUTPUT ARGUMENTS : */
7fd59977 2729/* ------------------- */
0d969553
Y
2730/* FPNTAB: Auxiliary table. */
2731/* VHERMT: Table of 2*(IORDRV+1) coeff. of 2*(IORDRV+1) polynoms of Hermite. */
2732/* SOSOTB: Table where the terms of constraints are added */
7fd59977 2733/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2734/* with ui and vj positive roots of the Legendre polynom */
2735/* of degree NBPNTU and NBPNTV respectively. */
2736/* DISOTB: Table where the terms of constraints are added */
7fd59977 2737/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2738/* with ui and vj positive roots of the polynom of Legendre */
2739/* of degree NBPNTU and NBPNTV respectively. */
2740/* SODITB: Table where the terms of constraints are added */
7fd59977 2741/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2742/* with ui and vj positive roots of the polynom of Legendre */
2743/* of degree NBPNTU and NBPNTV respectively. */
2744/* DIDITB: Table where the terms of constraints are added */
7fd59977 2745/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2746/* with ui and vj positive roots of the polynom of Legendre */
2747/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 2748
0d969553 2749/* COMMONS USED : */
7fd59977 2750/* ---------------- */
2751
0d969553 2752/* REFERENCES CALLED : */
7fd59977 2753/* ----------------------- */
2754
0d969553 2755/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2756/* ----------------------------------- */
2757
2758
7fd59977 2759/* > */
2760/* **********************************************************************
2761*/
2762
0d969553 2763/* Name of the routine */
7fd59977 2764
2765
2766 /* Parameter adjustments */
2767 diditb_dim1 = *nbpntu / 2 + 1;
2768 diditb_dim2 = *nbpntv / 2 + 1;
2769 diditb_offset = diditb_dim1 * diditb_dim2;
2770 diditb -= diditb_offset;
2771 disotb_dim1 = *nbpntu / 2;
2772 disotb_dim2 = *nbpntv / 2;
2773 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
2774 disotb -= disotb_offset;
2775 soditb_dim1 = *nbpntu / 2;
2776 soditb_dim2 = *nbpntv / 2;
2777 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
2778 soditb -= soditb_offset;
2779 sosotb_dim1 = *nbpntu / 2 + 1;
2780 sosotb_dim2 = *nbpntv / 2 + 1;
2781 sosotb_offset = sosotb_dim1 * sosotb_dim2;
2782 sosotb -= sosotb_offset;
2783 --vrootl;
2784 vhermt_dim1 = (*iordrv << 1) + 2;
2785 vhermt_offset = vhermt_dim1;
2786 vhermt -= vhermt_offset;
2787 fpntab_dim1 = *nbpntv;
2788 fpntab_offset = fpntab_dim1 + 1;
2789 fpntab -= fpntab_offset;
2790 ditbv2_dim1 = *nbpntu / 2 + 1;
2791 ditbv2_dim2 = *ndimen;
2792 ditbv2_offset = ditbv2_dim1 * (ditbv2_dim2 + 1);
2793 ditbv2 -= ditbv2_offset;
2794 ditbv1_dim1 = *nbpntu / 2 + 1;
2795 ditbv1_dim2 = *ndimen;
2796 ditbv1_offset = ditbv1_dim1 * (ditbv1_dim2 + 1);
2797 ditbv1 -= ditbv1_offset;
2798 sotbv2_dim1 = *nbpntu / 2 + 1;
2799 sotbv2_dim2 = *ndimen;
2800 sotbv2_offset = sotbv2_dim1 * (sotbv2_dim2 + 1);
2801 sotbv2 -= sotbv2_offset;
2802 sotbv1_dim1 = *nbpntu / 2 + 1;
2803 sotbv1_dim2 = *ndimen;
2804 sotbv1_offset = sotbv1_dim1 * (sotbv1_dim2 + 1);
2805 sotbv1 -= sotbv1_offset;
2806
2807 /* Function Body */
2808 ibb = AdvApp2Var_SysBase::mnfndeb_();
2809 if (ibb >= 3) {
2810 AdvApp2Var_SysBase::mgenmsg_("MMA2CD2", 7L);
2811 }
2812
0d969553 2813/* ------------------- Discretization of Hermit polynoms -----------
7fd59977 2814*/
2815
2816 ncfhv = (*iordrv + 1) << 1;
2817 i__1 = ncfhv;
2818 for (ii = 1; ii <= i__1; ++ii) {
2819 i__2 = *nbpntv;
2820 for (jj = 1; jj <= i__2; ++jj) {
2821 AdvApp2Var_MathBase::mmmpocur_(&ncfhv, &c__1, &ncfhv, &vhermt[ii * vhermt_dim1], &
2822 vrootl[jj], &fpntab[jj + ii * fpntab_dim1]);
2823/* L60: */
2824 }
2825/* L50: */
2826 }
2827
0d969553 2828/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 2829*/
2830
2831 nuroo = *nbpntu / 2;
2832 nvroo = *nbpntv / 2;
2833
2834 i__1 = *ndimen;
2835 for (nd = 1; nd <= i__1; ++nd) {
2836 i__2 = *iordrv + 1;
2837 for (ii = 1; ii <= i__2; ++ii) {
2838
2839 i__3 = nuroo;
2840 for (kk = 1; kk <= i__3; ++kk) {
2841 bid1 = sotbv1[kk + (nd + ii * sotbv1_dim2) * sotbv1_dim1];
2842 bid2 = sotbv2[kk + (nd + ii * sotbv2_dim2) * sotbv2_dim1];
2843 bid3 = ditbv1[kk + (nd + ii * ditbv1_dim2) * ditbv1_dim1];
2844 bid4 = ditbv2[kk + (nd + ii * ditbv2_dim2) * ditbv2_dim1];
2845 i__4 = nvroo;
2846 for (jj = 1; jj <= i__4; ++jj) {
2847 jjp = (*nbpntv + 1) / 2 + jj;
2848 jjm = nvroo - jj + 1;
2849 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1] =
2850 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1]
2851 - bid1 * (fpntab[jjp + ((ii << 1) - 1) *
2852 fpntab_dim1] + fpntab[jjm + ((ii << 1) - 1) *
2853 fpntab_dim1]) - bid2 * (fpntab[jjp + (ii << 1) *
2854 fpntab_dim1] + fpntab[jjm + (ii << 1) *
2855 fpntab_dim1]);
2856 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1] =
2857 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1]
2858 - bid3 * (fpntab[jjp + ((ii << 1) - 1) *
2859 fpntab_dim1] + fpntab[jjm + ((ii << 1) - 1) *
2860 fpntab_dim1]) - bid4 * (fpntab[jjp + (ii << 1) *
2861 fpntab_dim1] + fpntab[jjm + (ii << 1) *
2862 fpntab_dim1]);
2863 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1] =
2864 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1]
2865 - bid1 * (fpntab[jjp + ((ii << 1) - 1) *
2866 fpntab_dim1] - fpntab[jjm + ((ii << 1) - 1) *
2867 fpntab_dim1]) - bid2 * (fpntab[jjp + (ii << 1) *
2868 fpntab_dim1] - fpntab[jjm + (ii << 1) *
2869 fpntab_dim1]);
2870 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1] =
2871 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1]
2872 - bid3 * (fpntab[jjp + ((ii << 1) - 1) *
2873 fpntab_dim1] - fpntab[jjm + ((ii << 1) - 1) *
2874 fpntab_dim1]) - bid4 * (fpntab[jjp + (ii << 1) *
2875 fpntab_dim1] - fpntab[jjm + (ii << 1) *
2876 fpntab_dim1]);
2877/* L400: */
2878 }
2879/* L300: */
2880 }
2881/* L200: */
2882 }
2883
0d969553
Y
2884/* ------------ Case when the discretization is done only on the roots */
2885/* ---------- of Legendre polynom of uneven degree, 0 is root */
2886
7fd59977 2887
2888 if (*nbpntv % 2 == 1) {
2889 i__2 = *iordrv + 1;
2890 for (ii = 1; ii <= i__2; ++ii) {
2891 i__3 = nuroo;
2892 for (kk = 1; kk <= i__3; ++kk) {
2893 bid1 = sotbv1[kk + (nd + ii * sotbv1_dim2) * sotbv1_dim1]
2894 * fpntab[nvroo + 1 + ((ii << 1) - 1) *
2895 fpntab_dim1] + sotbv2[kk + (nd + ii * sotbv2_dim2)
2896 * sotbv2_dim1] * fpntab[nvroo + 1 + (ii << 1) *
2897 fpntab_dim1];
2898 sosotb[kk + nd * sosotb_dim2 * sosotb_dim1] -= bid1;
2899 bid2 = ditbv1[kk + (nd + ii * ditbv1_dim2) * ditbv1_dim1]
2900 * fpntab[nvroo + 1 + ((ii << 1) - 1) *
2901 fpntab_dim1] + ditbv2[kk + (nd + ii * ditbv2_dim2)
2902 * ditbv2_dim1] * fpntab[nvroo + 1 + (ii << 1) *
2903 fpntab_dim1];
2904 diditb[kk + nd * diditb_dim2 * diditb_dim1] -= bid2;
2905/* L550: */
2906 }
2907/* L500: */
2908 }
2909 }
2910
2911 if (*nbpntu % 2 == 1) {
2912 i__2 = *iordrv + 1;
2913 for (ii = 1; ii <= i__2; ++ii) {
2914 i__3 = nvroo;
2915 for (jj = 1; jj <= i__3; ++jj) {
2916 jjp = (*nbpntv + 1) / 2 + jj;
2917 jjm = nvroo - jj + 1;
2918 bid1 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * (
2919 fpntab[jjp + ((ii << 1) - 1) * fpntab_dim1] +
2920 fpntab[jjm + ((ii << 1) - 1) * fpntab_dim1]) +
2921 sotbv2[(nd + ii * sotbv2_dim2) * sotbv2_dim1] * (
2922 fpntab[jjp + (ii << 1) * fpntab_dim1] + fpntab[
2923 jjm + (ii << 1) * fpntab_dim1]);
2924 sosotb[(jj + nd * sosotb_dim2) * sosotb_dim1] -= bid1;
2925 bid2 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * (
2926 fpntab[jjp + ((ii << 1) - 1) * fpntab_dim1] -
2927 fpntab[jjm + ((ii << 1) - 1) * fpntab_dim1]) +
2928 sotbv2[(nd + ii * sotbv2_dim2) * sotbv2_dim1] * (
2929 fpntab[jjp + (ii << 1) * fpntab_dim1] - fpntab[
2930 jjm + (ii << 1) * fpntab_dim1]);
2931 diditb[jj + nd * diditb_dim2 * diditb_dim1] -= bid2;
2932/* L650: */
2933 }
2934/* L600: */
2935 }
2936 }
2937
2938 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
2939 i__2 = *iordrv + 1;
2940 for (ii = 1; ii <= i__2; ++ii) {
2941 bid1 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * fpntab[
2942 nvroo + 1 + ((ii << 1) - 1) * fpntab_dim1] + sotbv2[(
2943 nd + ii * sotbv2_dim2) * sotbv2_dim1] * fpntab[nvroo
2944 + 1 + (ii << 1) * fpntab_dim1];
2945 sosotb[nd * sosotb_dim2 * sosotb_dim1] -= bid1;
2946/* L700: */
2947 }
2948 }
2949
2950/* L100: */
2951 }
2952 goto L9999;
2953
2954/* ------------------------------ The End -------------------------------
2955*/
2956
2957L9999:
2958 if (ibb >= 3) {
2959 AdvApp2Var_SysBase::mgsomsg_("MMA2CD2", 7L);
2960 }
2961 return 0;
2962} /* mma2cd2_ */
2963
2964//=======================================================================
2965//function : mma2cd3_
2966//purpose :
2967//=======================================================================
2968int mma2cd3_(integer *ndimen,
2969 integer *nbpntu,
2970 doublereal *urootl,
2971 integer *nbpntv,
2972 integer *iordru,
2973 doublereal *sotbu1,
2974 doublereal *sotbu2,
2975 doublereal *ditbu1,
2976 doublereal *ditbu2,
2977 doublereal *fpntab,
2978 doublereal *uhermt,
2979 doublereal *sosotb,
2980 doublereal *soditb,
2981 doublereal *disotb,
2982 doublereal *diditb)
2983
2984{
1ef32e96 2985 integer c__1 = 1;
41194117 2986
7fd59977 2987 /* System generated locals */
2988 integer sotbu1_dim1, sotbu1_dim2, sotbu1_offset, sotbu2_dim1, sotbu2_dim2,
2989 sotbu2_offset, ditbu1_dim1, ditbu1_dim2, ditbu1_offset,
2990 ditbu2_dim1, ditbu2_dim2, ditbu2_offset, fpntab_dim1,
2991 fpntab_offset, uhermt_dim1, uhermt_offset, sosotb_dim1,
2992 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2993 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2994 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4;
2995
2996 /* Local variables */
1ef32e96
RL
2997 integer ncfhu, nuroo, nvroo, ii, nd, jj, kk, ibb, kkm, kkp;
2998 doublereal bid1, bid2, bid3, bid4;
7fd59977 2999
3000/* **********************************************************************
3001*/
0d969553 3002/* FUNCTION : */
7fd59977 3003/* ---------- */
0d969553
Y
3004/* Discretisation on the parameters of polynoms of interpolation */
3005/* of constraints on 2 borders iso-U of order IORDRU. */
7fd59977 3006
0d969553
Y
3007
3008/* KEYWORDS : */
7fd59977 3009/* ----------- */
3010/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
3011
0d969553 3012/* INPUT ARGUMENTS : */
7fd59977 3013/* ------------------ */
0d969553
Y
3014/* NDIMEN: Dimension of the space. */
3015/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
3016/* This is also the nb of root of Legendre polynom where discretization is done. */
3017/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
3018*/
3019/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
3020/* This is also the nb of root of Legendre polynom where discretization is done. */
3021/* IORDRV: Order of constraint imposed at the extremities of iso-V */
3022/* = 0, calculate the extremities of iso-V */
3023/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
3024/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
3025/* SOTBU1: Table of NBPNTU/2 sums of 2 index points */
3026/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V0. */
3027/* SOTBU2: Table of NBPNTV/2 sums of 2 index points */
3028/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V1. */
3029/* DITBU1: Table of NBPNTU/2 differences of 2 index points */
3030/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V0. */
3031/* DITBU2: Table of NBPNTU/2 differences of 2 index points */
3032/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V1. */
3033/* SOSOTB: Preinitialized table (input/output argument). */
3034/* DISOTB: Preinitialized table (input/output argument). */
3035/* SODITB: Preinitialized table (input/output argument). */
3036/* DIDITB: Preinitialized table (input/output argument) */
3037
3038/* OUTPUT ARGUMENTS : */
7fd59977 3039/* ------------------- */
0d969553
Y
3040/* FPNTAB: Auxiliary table. */
3041/* UHERMT: Table of 2*(IORDRU+1) coeff. of 2*(IORDRU+1) polynoms of Hermite. */
3042/* SOSOTB: Table where the terms of constraints are added */
7fd59977 3043/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
3044/* with ui and vj positive roots of the Legendre polynom */
3045/* of degree NBPNTU and NBPNTV respectively. */
3046/* DISOTB: Table where the terms of constraints are added */
7fd59977 3047/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
3048/* with ui and vj positive roots of the polynom of Legendre */
3049/* of degree NBPNTU and NBPNTV respectively. */
3050/* SODITB: Table where the terms of constraints are added */
7fd59977 3051/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
3052/* with ui and vj positive roots of the polynom of Legendre */
3053/* of degree NBPNTU and NBPNTV respectively. */
3054/* DIDITB: Table where the terms of constraints are added */
7fd59977 3055/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
3056/* with ui and vj positive roots of the polynom of Legendre */
3057/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 3058
0d969553 3059/* COMMONS USED : */
7fd59977 3060/* ---------------- */
3061
0d969553 3062/* REFERENCES CALLED : */
7fd59977 3063/* ----------------------- */
3064
0d969553 3065/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 3066/* ----------------------------------- */
3067
7fd59977 3068/* $ HISTORIQUE DES MODIFICATIONS : */
3069/* -------------------------------- */
3070/* 08-08-1991: RBD; Creation. */
3071/* > */
3072/* **********************************************************************
3073*/
3074
0d969553 3075/* Name of the routine */
7fd59977 3076
3077
3078 /* Parameter adjustments */
3079 --urootl;
3080 diditb_dim1 = *nbpntu / 2 + 1;
3081 diditb_dim2 = *nbpntv / 2 + 1;
3082 diditb_offset = diditb_dim1 * diditb_dim2;
3083 diditb -= diditb_offset;
3084 disotb_dim1 = *nbpntu / 2;
3085 disotb_dim2 = *nbpntv / 2;
3086 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
3087 disotb -= disotb_offset;
3088 soditb_dim1 = *nbpntu / 2;
3089 soditb_dim2 = *nbpntv / 2;
3090 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
3091 soditb -= soditb_offset;
3092 sosotb_dim1 = *nbpntu / 2 + 1;
3093 sosotb_dim2 = *nbpntv / 2 + 1;
3094 sosotb_offset = sosotb_dim1 * sosotb_dim2;
3095 sosotb -= sosotb_offset;
3096 uhermt_dim1 = (*iordru << 1) + 2;
3097 uhermt_offset = uhermt_dim1;
3098 uhermt -= uhermt_offset;
3099 fpntab_dim1 = *nbpntu;
3100 fpntab_offset = fpntab_dim1 + 1;
3101 fpntab -= fpntab_offset;
3102 ditbu2_dim1 = *nbpntv / 2 + 1;
3103 ditbu2_dim2 = *ndimen;
3104 ditbu2_offset = ditbu2_dim1 * (ditbu2_dim2 + 1);
3105 ditbu2 -= ditbu2_offset;
3106 ditbu1_dim1 = *nbpntv / 2 + 1;
3107 ditbu1_dim2 = *ndimen;
3108 ditbu1_offset = ditbu1_dim1 * (ditbu1_dim2 + 1);
3109 ditbu1 -= ditbu1_offset;
3110 sotbu2_dim1 = *nbpntv / 2 + 1;
3111 sotbu2_dim2 = *ndimen;
3112 sotbu2_offset = sotbu2_dim1 * (sotbu2_dim2 + 1);
3113 sotbu2 -= sotbu2_offset;
3114 sotbu1_dim1 = *nbpntv / 2 + 1;
3115 sotbu1_dim2 = *ndimen;
3116 sotbu1_offset = sotbu1_dim1 * (sotbu1_dim2 + 1);
3117 sotbu1 -= sotbu1_offset;
3118
3119 /* Function Body */
3120 ibb = AdvApp2Var_SysBase::mnfndeb_();
3121 if (ibb >= 3) {
3122 AdvApp2Var_SysBase::mgenmsg_("MMA2CD3", 7L);
3123 }
3124
0d969553 3125/* ------------------- Discretization of polynoms of Hermit -----------
7fd59977 3126*/
3127
3128 ncfhu = (*iordru + 1) << 1;
3129 i__1 = ncfhu;
3130 for (ii = 1; ii <= i__1; ++ii) {
3131 i__2 = *nbpntu;
3132 for (kk = 1; kk <= i__2; ++kk) {
3133 AdvApp2Var_MathBase::mmmpocur_(&ncfhu,
3134 &c__1,
3135 &ncfhu,
3136 &uhermt[ii * uhermt_dim1],
3137 &urootl[kk],
3138 &fpntab[kk + ii * fpntab_dim1]);
3139/* L60: */
3140 }
3141/* L50: */
3142 }
3143
0d969553 3144/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 3145*/
3146
3147 nvroo = *nbpntv / 2;
3148 nuroo = *nbpntu / 2;
3149
3150 i__1 = *ndimen;
3151 for (nd = 1; nd <= i__1; ++nd) {
3152 i__2 = *iordru + 1;
3153 for (ii = 1; ii <= i__2; ++ii) {
3154
3155 i__3 = nvroo;
3156 for (jj = 1; jj <= i__3; ++jj) {
3157 bid1 = sotbu1[jj + (nd + ii * sotbu1_dim2) * sotbu1_dim1];
3158 bid2 = sotbu2[jj + (nd + ii * sotbu2_dim2) * sotbu2_dim1];
3159 bid3 = ditbu1[jj + (nd + ii * ditbu1_dim2) * ditbu1_dim1];
3160 bid4 = ditbu2[jj + (nd + ii * ditbu2_dim2) * ditbu2_dim1];
3161 i__4 = nuroo;
3162 for (kk = 1; kk <= i__4; ++kk) {
3163 kkp = (*nbpntu + 1) / 2 + kk;
3164 kkm = nuroo - kk + 1;
3165 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1] =
3166 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1]
3167 - bid1 * (fpntab[kkp + ((ii << 1) - 1) *
3168 fpntab_dim1] + fpntab[kkm + ((ii << 1) - 1) *
3169 fpntab_dim1]) - bid2 * (fpntab[kkp + (ii << 1) *
3170 fpntab_dim1] + fpntab[kkm + (ii << 1) *
3171 fpntab_dim1]);
3172 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1] =
3173 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1]
3174 - bid1 * (fpntab[kkp + ((ii << 1) - 1) *
3175 fpntab_dim1] - fpntab[kkm + ((ii << 1) - 1) *
3176 fpntab_dim1]) - bid2 * (fpntab[kkp + (ii << 1) *
3177 fpntab_dim1] - fpntab[kkm + (ii << 1) *
3178 fpntab_dim1]);
3179 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1] =
3180 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1]
3181 - bid3 * (fpntab[kkp + ((ii << 1) - 1) *
3182 fpntab_dim1] + fpntab[kkm + ((ii << 1) - 1) *
3183 fpntab_dim1]) - bid4 * (fpntab[kkp + (ii << 1) *
3184 fpntab_dim1] + fpntab[kkm + (ii << 1) *
3185 fpntab_dim1]);
3186 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1] =
3187 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1]
3188 - bid3 * (fpntab[kkp + ((ii << 1) - 1) *
3189 fpntab_dim1] - fpntab[kkm + ((ii << 1) - 1) *
3190 fpntab_dim1]) - bid4 * (fpntab[kkp + (ii << 1) *
3191 fpntab_dim1] - fpntab[kkm + (ii << 1) *
3192 fpntab_dim1]);
3193/* L400: */
3194 }
3195/* L300: */
3196 }
3197/* L200: */
3198 }
3199
0d969553
Y
3200/* ------------ Case when the discretization is done only on the roots */
3201/* ---------- of Legendre polynom of uneven degree, 0 is root */
3202
3203
7fd59977 3204
3205 if (*nbpntu % 2 == 1) {
3206 i__2 = *iordru + 1;
3207 for (ii = 1; ii <= i__2; ++ii) {
3208 i__3 = nvroo;
3209 for (jj = 1; jj <= i__3; ++jj) {
3210 bid1 = sotbu1[jj + (nd + ii * sotbu1_dim2) * sotbu1_dim1]
3211 * fpntab[nuroo + 1 + ((ii << 1) - 1) *
3212 fpntab_dim1] + sotbu2[jj + (nd + ii * sotbu2_dim2)
3213 * sotbu2_dim1] * fpntab[nuroo + 1 + (ii << 1) *
3214 fpntab_dim1];
3215 sosotb[(jj + nd * sosotb_dim2) * sosotb_dim1] -= bid1;
3216 bid2 = ditbu1[jj + (nd + ii * ditbu1_dim2) * ditbu1_dim1]
3217 * fpntab[nuroo + 1 + ((ii << 1) - 1) *
3218 fpntab_dim1] + ditbu2[jj + (nd + ii * ditbu2_dim2)
3219 * ditbu2_dim1] * fpntab[nuroo + 1 + (ii << 1) *
3220 fpntab_dim1];
3221 diditb[(jj + nd * diditb_dim2) * diditb_dim1] -= bid2;
3222/* L550: */
3223 }
3224/* L500: */
3225 }
3226 }
3227
3228 if (*nbpntv % 2 == 1) {
3229 i__2 = *iordru + 1;
3230 for (ii = 1; ii <= i__2; ++ii) {
3231 i__3 = nuroo;
3232 for (kk = 1; kk <= i__3; ++kk) {
3233 kkp = (*nbpntu + 1) / 2 + kk;
3234 kkm = nuroo - kk + 1;
3235 bid1 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * (
3236 fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] +
3237 fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1]) +
3238 sotbu2[(nd + ii * sotbu2_dim2) * sotbu2_dim1] * (
3239 fpntab[kkp + (ii << 1) * fpntab_dim1] + fpntab[
3240 kkm + (ii << 1) * fpntab_dim1]);
3241 sosotb[kk + nd * sosotb_dim2 * sosotb_dim1] -= bid1;
3242 bid2 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * (
3243 fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] -
3244 fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1]) +
3245 sotbu2[(nd + ii * sotbu2_dim2) * sotbu2_dim1] * (
3246 fpntab[kkp + (ii << 1) * fpntab_dim1] - fpntab[
3247 kkm + (ii << 1) * fpntab_dim1]);
3248 diditb[kk + nd * diditb_dim2 * diditb_dim1] -= bid2;
3249/* L650: */
3250 }
3251/* L600: */
3252 }
3253 }
3254
3255 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
3256 i__2 = *iordru + 1;
3257 for (ii = 1; ii <= i__2; ++ii) {
3258 bid1 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * fpntab[
3259 nuroo + 1 + ((ii << 1) - 1) * fpntab_dim1] + sotbu2[(
3260 nd + ii * sotbu2_dim2) * sotbu2_dim1] * fpntab[nuroo
3261 + 1 + (ii << 1) * fpntab_dim1];
3262 sosotb[nd * sosotb_dim2 * sosotb_dim1] -= bid1;
3263/* L700: */
3264 }
3265 }
3266
3267/* L100: */
3268 }
3269 goto L9999;
3270
3271/* ------------------------------ The End -------------------------------
3272*/
3273
3274L9999:
3275 if (ibb >= 3) {
3276 AdvApp2Var_SysBase::mgsomsg_("MMA2CD3", 7L);
3277 }
3278 return 0;
3279} /* mma2cd3_ */
3280
3281//=======================================================================
3282//function : mma2cdi_
3283//purpose :
3284//=======================================================================
3285int AdvApp2Var_ApproxF2var::mma2cdi_( integer *ndimen,
3286 integer *nbpntu,
3287 doublereal *urootl,
3288 integer *nbpntv,
3289 doublereal *vrootl,
3290 integer *iordru,
3291 integer *iordrv,
3292 doublereal *contr1,
3293 doublereal *contr2,
3294 doublereal *contr3,
3295 doublereal *contr4,
3296 doublereal *sotbu1,
3297 doublereal *sotbu2,
3298 doublereal *ditbu1,
3299 doublereal *ditbu2,
3300 doublereal *sotbv1,
3301 doublereal *sotbv2,
3302 doublereal *ditbv1,
3303 doublereal *ditbv2,
3304 doublereal *sosotb,
3305 doublereal *soditb,
3306 doublereal *disotb,
3307 doublereal *diditb,
3308 integer *iercod)
3309
3310{
1ef32e96 3311 integer c__8 = 8;
7fd59977 3312
3313 /* System generated locals */
3314 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
3315 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
3316 contr4_dim1, contr4_dim2, contr4_offset, sosotb_dim1, sosotb_dim2,
3317 sosotb_offset, diditb_dim1, diditb_dim2, diditb_offset,
3318 soditb_dim1, soditb_dim2, soditb_offset, disotb_dim1, disotb_dim2,
3319 disotb_offset;
3320
3321 /* Local variables */
1ef32e96
RL
3322 integer ilong;
3323 intptr_t iofwr;
3324 doublereal* wrkar = 0;