1 /* Siconos is a program dedicated to modeling, simulation and control
2  * of non smooth dynamical systems.
3  *
4  * Copyright 2021 INRIA.
5  *
6  * Licensed under the Apache License, Version 2.0 (the "License");
7  * you may not use this file except in compliance with the License.
8  * You may obtain a copy of the License at
9  *
10  * http://www.apache.org/licenses/LICENSE-2.0
11  *
12  * Unless required by applicable law or agreed to in writing, software
13  * distributed under the License is distributed on an "AS IS" BASIS,
14  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15  * See the License for the specific language governing permissions and
16  * limitations under the License.
17 */
18 
19 #include <stdio.h>                      // for NULL, fprintf, stderr
20 #include <stdlib.h>                     // for exit, EXIT_FAILURE
21 #include "Friction_cst.h"               // for SICONOS_FRICTION_3D_NCPGlocke...
22 #include "NCP_FixedP.h"                 // for Fixe
23 #include "NumericsFwd.h"                // for SolverOptions, FrictionContac...
24 #include "SolverOptions.h"              // for SolverOptions
25 #include "fc3d_2NCP_Glocker.h"          // for NCPGlocker_initialize, comput...
26 #include "fc3d_NCPGlockerFixedPoint.h"  // for F_GlockerFixedP, fc3d_FixedP_...
27 #include "fc3d_Solvers.h"               // for FreeSolverPtr, PostSolverPtr
28 #include "SiconosBlas.h"                      // for cblas_dcopy
29 
30 /* Pointer to function used to update the solver, to formalize the local problem for example. */
31 typedef void (*UpdateSolverPtr)(int, double*);
32 
33 static UpdateSolverPtr updateSolver = NULL;
34 static PostSolverPtr postSolver = NULL;
35 static FreeSolverPtr freeSolver = NULL;
36 
37 /* size of a block */
38 static int Fsize;
39 
40 /** writes \f$ F(z) \f$ using Glocker formulation
41  */
F_GlockerFixedP(int sizeF,double * reaction,double * FVector,int up2Date)42 void F_GlockerFixedP(int sizeF, double* reaction, double* FVector, int up2Date)
43 {
44   /* Glocker formulation */
45   double* FGlocker = NULL;
46   computeFGlocker(&FGlocker, up2Date);
47   /* Note that FGlocker is a static var. in fc3d2NCP_Glocker and thus there is no memory allocation in
48      the present file.
49   */
50 
51   /* TMP COPY: review memory management for FGlocker ...*/
52   cblas_dcopy(sizeF, FGlocker, 1, FVector, 1);
53   FGlocker = NULL;
54 }
55 
56 /** writes \f$ \nabla_z F(z) \f$  using Glocker formulation and the Fischer-Burmeister function.
57  */
58 
59 
fc3d_FixedP_initialize(FrictionContactProblem * problem,FrictionContactProblem * localproblem,SolverOptions * localsolver_options)60 void fc3d_FixedP_initialize(FrictionContactProblem* problem, FrictionContactProblem* localproblem, SolverOptions * localsolver_options)
61 {
62 
63   /*
64      Initialize solver (Compute F) according to the chosen formulation.
65   */
66 
67   /* Glocker formulation */
68   if(localsolver_options->solverId == SICONOS_FRICTION_3D_NCPGlockerFBFixedPoint)
69   {
70     Fsize = 5;
71     NCPGlocker_initialize(problem, localproblem);
72     /*     updateSolver = &NCPGlocker_update; */
73     postSolver = &NCPGlocker_post;
74     freeSolver = &NCPGlocker_free;
75   }
76   else
77   {
78     fprintf(stderr, "Numerics, fc3d_nsgs failed. Unknown formulation type.\n");
79     exit(EXIT_FAILURE);
80   }
81 }
82 
fc3d_FixedP_solve(FrictionContactProblem * localproblem,double * reaction,SolverOptions * options)83 int fc3d_FixedP_solve(FrictionContactProblem * localproblem, double* reaction, SolverOptions * options)
84 {
85   int * iparam = options->iparam;
86   double * dparam = options->dparam;
87 
88   double * reactionBlock = reaction;
89 
90   int info = Fixe(Fsize, reactionBlock, iparam, dparam);
91 
92   if(info > 0)
93   {
94     fprintf(stderr, "Numerics, fc3d_FixedP failed, reached max. number of iterations without convergence. Residual = %f\n", dparam[SICONOS_DPARAM_RESIDU]);
95     exit(EXIT_FAILURE);
96   }
97   return info;
98 
99   /*   (*postSolver)(contact,reaction); */
100 }
101 
fc3d_FixedP_free(FrictionContactProblem * problem,FrictionContactProblem * localproblem,SolverOptions * localsolver_option)102 void fc3d_FixedP_free(FrictionContactProblem * problem, FrictionContactProblem * localproblem, SolverOptions * localsolver_option)
103 {
104   updateSolver = NULL;
105   postSolver = NULL;
106   (*freeSolver)();
107 }
108 
109 /*
110 double fc3d_FixedP_computeError(int contact, int dimReaction, double* reaction, double * error)
111 {
112   return 0.0;
113 }
114 */
115