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 #include "SiconosConfig.h" // for WITH_FCLIB // IWYU pragma: keep
19
20 #ifdef WITH_FCLIB
21 #define DEBUG_NOCOLOR
22 #define DEBUG_STDOUT
23 #define DEBUG_MESSAGES
24 #include "fclib_interface.h"
25 #include <assert.h> // for assert
26 #include <fclib.h> // for fclib_matrix, fclib_global
27 #include <stdio.h> // for NULL, fprintf, stderr
28 #include <stdlib.h> // for malloc, free, exit, EXIT_F...
29 #include "CSparseMatrix_internal.h" // for CSparseMatrix, CS_INT, cs_...
30 #include "FrictionContactProblem.h" // for FrictionContactProblem
31 #include "GlobalFrictionContactProblem.h" // for GlobalFrictionContactProblem
32 #include "GlobalRollingFrictionContactProblem.h" // for GlobalRollingFrictionContactProblem
33 #include "NumericsMatrix.h" // for NumericsMatrix, RawNumeric...
34 #include "NumericsSparseMatrix.h" // for NumericsSparseMatrix, NSM_CSC
35 #include "SiconosConfig.h" // for WITH_FCLIB
36 #include "SparseBlockMatrix.h" // for SBM_from_csparse, SBM_to_s...
37 #include "siconos_debug.h" // for DEBUG_PRINT, DEBUG_PRINTF
38 #include "timers_interf.h" // for MAYBE_UNUSED
39
40 // avoid a conflict with old csparse.h in case fclib includes it
41 #define _CS_H
42
43
int_to_csi(int * o,CS_INT * d,unsigned int n)44 static void int_to_csi(int* o, CS_INT* d, unsigned int n)
45 {
46 for(unsigned int i=0; i<n; ++i)
47 {
48 d[i] = (CS_INT) o[i];
49 }
50 }
51
csi_to_int(CS_INT * o,int * d,unsigned int n)52 static void csi_to_int(CS_INT* o, int* d, unsigned int n)
53 {
54 for(unsigned int i=0; i<n; ++i)
55 {
56 d[i] = (int) o[i];
57 }
58 }
59
60
61
from_fclib_local(const fclib_local * fclib_problem)62 FrictionContactProblem* from_fclib_local(const fclib_local* fclib_problem)
63 {
64 FrictionContactProblem* problem;
65
66 problem = (FrictionContactProblem*)malloc(sizeof(FrictionContactProblem));
67
68 problem->dimension = fclib_problem->spacedim;
69 problem->mu = fclib_problem->mu;
70 problem->q = fclib_problem->q;
71
72 problem->numberOfContacts = fclib_problem->W->m / fclib_problem->spacedim; /* cf fclib spec */
73
74 problem->M = NM_create(NM_SPARSE_BLOCK, fclib_problem->W->m, fclib_problem->W->n);
75
76 problem->M->matrix1->block = NULL;
77 problem->M->matrix1->index1_data = NULL;
78 problem->M->matrix1->index2_data = NULL;
79
80 CSparseMatrix W;
81
82 W.nzmax = (CS_INT) fclib_problem->W->nzmax;
83 W.m = (CS_INT) fclib_problem->W->m;
84 W.n = (CS_INT) fclib_problem->W->n;
85
86 if(fclib_problem->W->nz == -1)
87 {
88 /* compressed colums */
89 W.p = (CS_INT*) malloc(sizeof(CS_INT)*(W.n+1));
90 int_to_csi(fclib_problem->W->p, W.p, (unsigned)(W.n+1));
91 }
92 else if(fclib_problem->W->nz == -2)
93 {
94 /* compressed rows */
95 W.p = (CS_INT*) malloc(sizeof(CS_INT)*(W.m+1));
96 int_to_csi(fclib_problem->W->p, W.p, (unsigned)(W.m+1));
97 }
98 else
99 {
100 /* triplet */
101 W.p = (CS_INT*) malloc(sizeof(CS_INT)*W.nzmax);
102 int_to_csi(fclib_problem->W->p, W.p, (unsigned) W.nzmax);
103 }
104
105 W.i = (CS_INT*) malloc(sizeof(CS_INT)*W.nzmax);
106 int_to_csi(fclib_problem->W->i, W.i, (unsigned) W.nzmax);
107
108 W.x = fclib_problem->W->x;
109
110 W.nz = fclib_problem->W->nz;
111
112 SBM_from_csparse(problem->dimension, &W, problem->M->matrix1);
113
114 free(W.p);
115 free(W.i);
116
117 NM_reset_versions(problem->M);
118
119 return problem;
120
121 }
122
123
frictionContact_fclib_read(const char * path)124 FrictionContactProblem* frictionContact_fclib_read(const char *path)
125 {
126
127 fclib_local *fclib_problem;
128
129 fclib_problem = fclib_read_local(path);
130
131 if(!fclib_problem)
132 {
133 return NULL;
134 }
135
136 return from_fclib_local(fclib_problem);
137 }
138
frictionContact_fclib_write(FrictionContactProblem * problem,char * title,char * description,char * mathInfo,const char * path,int ndof)139 int frictionContact_fclib_write(FrictionContactProblem* problem, char * title, char * description, char * mathInfo,
140 const char *path, int ndof)
141 {
142 int info = 0;
143
144 fclib_local *fclib_problem;
145
146 fclib_problem = (fclib_local*)malloc(sizeof(fclib_local));
147
148 fclib_problem->spacedim = problem->dimension;
149 fclib_problem->mu = problem->mu;
150 fclib_problem->q = problem->q;
151
152 fclib_problem->s = NULL;
153
154 fclib_problem->info = (struct fclib_info*)malloc(sizeof(struct fclib_info)) ;
155 fclib_problem->info->title = title;
156 fclib_problem->info->description = description;
157 fclib_problem->info->math_info = mathInfo;
158
159 fclib_problem->W = (struct fclib_matrix*)malloc(sizeof(struct fclib_matrix));
160 fclib_problem->R = NULL;
161 fclib_problem->V = NULL;
162
163 fclib_problem->W->m = problem->M->size0;
164 fclib_problem->W->n = problem->M->size1;
165
166
167 CSparseMatrix * spmat = NULL;
168
169 if(problem ->M->storageType == NM_DENSE) /* Dense Matrix */
170 {
171 fclib_problem->W->nzmax = problem->M->size0 * problem->M->size1;
172 fclib_problem->W->p = (int*)malloc((fclib_problem->W->m + 1) * sizeof(int));
173 fclib_problem->W->i = (int*)malloc((fclib_problem->W->nzmax) * sizeof(int));
174 fclib_problem->W->x = (double*)malloc((fclib_problem->W->nzmax) * sizeof(double));
175 fclib_problem->W->nz = -2;
176 fclib_problem->W->info = NULL;
177 for(int i = 0; i < problem ->M->size0 ; i++)
178 {
179 fclib_problem->W->p[i] = i * problem ->M->size1;
180 for(int j = 0; j < problem ->M->size1 ; j++)
181 {
182 fclib_problem->W->x[i * problem ->M->size1 + j ] = problem ->M->matrix0[j * problem ->M->size0 + i ];
183 fclib_problem->W->i[i * problem ->M->size1 + j ] = j;
184 }
185 }
186 fclib_problem->W->p[fclib_problem->W->m] = (fclib_problem->W->m) * problem ->M->size1;
187
188 }
189 else if(problem ->M->storageType == NM_SPARSE_BLOCK) /* Sparse block storage */
190 {
191 spmat = (CSparseMatrix*)malloc(sizeof(CSparseMatrix));
192 int MAYBE_UNUSED res = SBM_to_sparse_init_memory(problem ->M->matrix1, spmat);
193 res = SBM_to_sparse(problem->M->matrix1, spmat);
194 fclib_problem->W->nzmax = (int) spmat->nzmax;
195 fclib_problem->W->m = (int) spmat->m;
196 fclib_problem->W->n = (int) spmat->n;
197 fclib_problem->W->x = spmat->x;
198 fclib_problem->W->nz = (int) spmat->nz;
199
200 if(spmat->nz == -1)
201 {
202 fclib_problem->W->p = (int*) malloc(sizeof(int)*(spmat->n+1));
203 csi_to_int(spmat->p, fclib_problem->W->p, (unsigned)(spmat->n+1));
204 }
205 else if(spmat->nz == -2)
206 {
207 fclib_problem->W->p = (int*) malloc(sizeof(int)*(spmat->m+1));
208 csi_to_int(spmat->p, fclib_problem->W->p, (unsigned)(spmat->m+1));
209 }
210 else
211 {
212 fclib_problem->W->p = (int*) malloc(sizeof(int)*spmat->nzmax);
213 csi_to_int(spmat->p, fclib_problem->W->p, (unsigned) spmat->nzmax);
214 }
215
216 fclib_problem->W->i = (int*) malloc(sizeof(int)*spmat->nzmax);
217 csi_to_int(spmat->i, fclib_problem->W->i, (unsigned) spmat->nzmax);
218
219 fclib_problem->W->info = NULL;
220 }
221 else
222 {
223 fprintf(stderr, "frictionContact_fclib_write, unknown storage type for A.\n");
224 exit(EXIT_FAILURE); ;
225 }
226
227 info = fclib_write_local(fclib_problem, path);
228
229 info = fclib_create_int_attributes_in_info(path, "numberOfDegreeOfFreedom",
230 ndof);
231
232
233
234 /* fclib_delete_local (fclib_problem); */
235
236 if(problem ->M->storageType == NM_DENSE) /* Dense Matrix */
237 {
238 free(fclib_problem->W->x);
239 }
240 else if(problem ->M->storageType == NM_SPARSE_BLOCK)
241 {
242 cs_spfree(spmat);
243 }
244 free(fclib_problem->W->p);
245 free(fclib_problem->W->i);
246 free(fclib_problem->W);
247 free(fclib_problem->info);
248
249 free(fclib_problem);
250
251 return info;
252
253 }
frictionContact_fclib_write_guess(double * reaction,double * velocity,const char * path)254 int frictionContact_fclib_write_guess(double * reaction, double * velocity,
255 const char *path)
256 {
257 int info = 0;
258 int number_of_guesses = 1;
259 fclib_solution *guesses = (fclib_solution *) malloc(number_of_guesses*sizeof(fclib_solution));
260 guesses->v = NULL;
261 guesses->l = NULL;
262 guesses->u = velocity;
263 guesses->r = reaction;
264
265 info = fclib_write_guesses(number_of_guesses, guesses, path);
266 return info;
267
268 }
269
from_fclib_global(const fclib_global * fclib_problem)270 GlobalFrictionContactProblem* from_fclib_global(const fclib_global* fclib_problem)
271 {
272 GlobalFrictionContactProblem* problem;
273
274 problem = (GlobalFrictionContactProblem*)malloc(sizeof(GlobalFrictionContactProblem));
275
276 problem->dimension = fclib_problem->spacedim;
277 problem->mu = fclib_problem->mu;
278 problem->q = fclib_problem->f;
279 problem->b = fclib_problem->w;
280 problem->env = NULL;
281
282 problem->numberOfContacts = fclib_problem->H->n / fclib_problem->spacedim; /* cf fclib spec */
283
284 problem->M = NM_create(NM_SPARSE, fclib_problem->M->m, fclib_problem->M->n);
285
286 CSparseMatrix * M = (CSparseMatrix*)malloc(sizeof(CSparseMatrix));
287 M->nzmax = (CS_INT) fclib_problem->M->nzmax;
288 M->m = (CS_INT) fclib_problem->M->m;
289 M->n = (CS_INT) fclib_problem->M->n;
290
291 M->x = fclib_problem->M->x;
292
293 if(fclib_problem->M->nz == -1)
294 {
295 /* compressed colums */
296 problem->M->matrix2->csc= M;
297 problem->M->matrix2->origin = NSM_CSC;
298 M->nz = (CS_INT) fclib_problem->M->nz;
299 M->p = (CS_INT*) malloc(sizeof(CS_INT)*(M->n+1));
300 int_to_csi(fclib_problem->M->p, M->p, (unsigned)(M->n+1));
301 }
302 else if(fclib_problem->M->nz == -2)
303 {
304 /* compressed rows */
305 M->nz = (CS_INT) fclib_problem->M->nz;
306 M->p = (CS_INT*) malloc(sizeof(CS_INT)*(M->m+1));
307 int_to_csi(fclib_problem->M->p, M->p, (unsigned)(M->m+1));
308 /* since problem->M->matrix2->csr does not exist, we need
309 to fill transform M into a triplet or csc before returning
310 */
311
312 fprintf(stderr, "from_fclib_local not implemented for csr matrices.\n");
313 exit(EXIT_FAILURE); ;
314 }
315 else
316 {
317 /* triplet */
318 problem->M->matrix2->triplet=M;
319 problem->M->matrix2->origin = NSM_TRIPLET;
320 M->nz = (CS_INT) fclib_problem->M->nz;
321 M->p = (CS_INT*) malloc(sizeof(CS_INT)*M->nzmax);
322 int_to_csi(fclib_problem->M->p, M->p, (unsigned) M->nz);
323 }
324 M->i = (CS_INT*) malloc(sizeof(CS_INT)*M->nzmax);
325 int_to_csi(fclib_problem->M->i, M->i, (unsigned) M->nz);
326
327
328 problem->H = NM_create(NM_SPARSE, fclib_problem->H->m, fclib_problem->H->n);
329
330 CSparseMatrix * H = (CSparseMatrix*)malloc(sizeof(CSparseMatrix));
331
332 H->nzmax = (CS_INT) fclib_problem->H->nzmax;
333 H->m = (CS_INT) fclib_problem->H->m;
334 H->n = (CS_INT) fclib_problem->H->n;
335 H->nz = (CS_INT) fclib_problem->H->nz;
336 H->x = fclib_problem->H->x;
337
338 if(fclib_problem->H->nz == -1)
339 {
340 /* compressed colums */
341 problem->H->matrix2->csc= H;
342 problem->H->matrix2->origin = NSM_CSC;
343 H->p = (CS_INT*) malloc(sizeof(CS_INT)*(H->n+1));
344 int_to_csi(fclib_problem->H->p, H->p, (unsigned)(H->n+1));
345 }
346 else if(fclib_problem->H->nz == -2)
347 {
348 /* compressed rows */
349 fprintf(stderr, "from_fclib_local not implemented for csr matrices.\n");
350 exit(EXIT_FAILURE); ;
351 }
352 else
353 {
354 /* triplet */
355 problem->H->matrix2->triplet=H;
356 problem->H->matrix2->origin = NSM_TRIPLET;
357 H->p = (CS_INT*) malloc(sizeof(CS_INT)*H->nzmax);
358 int_to_csi(fclib_problem->H->p, H->p, (unsigned) H->nz);
359 }
360
361 H->i = (CS_INT*) malloc(sizeof(CS_INT)*H->nzmax);
362 int_to_csi(fclib_problem->H->i, H->i, (unsigned) H->nz);
363
364
365 NM_reset_versions(problem->M);
366 return problem;
367
368 }
369
370
globalFrictionContact_fclib_read(const char * path)371 GlobalFrictionContactProblem* globalFrictionContact_fclib_read(const char *path)
372 {
373
374 fclib_global *fclib_problem;
375
376 fclib_problem = fclib_read_global(path);
377
378 if(!fclib_problem)
379 {
380 return NULL;
381 }
382
383 return from_fclib_global(fclib_problem);
384 }
385
386
globalFrictionContact_fclib_write(GlobalFrictionContactProblem * problem,char * title,char * description,char * mathInfo,const char * path)387 int globalFrictionContact_fclib_write(
388 GlobalFrictionContactProblem* problem,
389 char * title,
390 char * description,
391 char * mathInfo,
392 const char *path)
393 {
394 int rinfo = 0;
395
396 DEBUG_PRINTF("construction of fclib_problem in %s with title = %s and description = %s\n", path, title, description);
397 if(problem->numberOfContacts == 0)
398 {
399 DEBUG_PRINT("zero contacts");
400 return rinfo;
401
402 }
403 /* globalFrictionContact_display(problem); */
404 /* FILE * file = fopen("toto.dat", "w"); */
405 /* globalFrictionContact_printInFile(problem, file); */
406
407 fclib_global *fclib_problem;
408 fclib_problem = (fclib_global*)malloc(sizeof(fclib_global));
409
410 fclib_problem->info = (struct fclib_info*)malloc(sizeof(struct fclib_info)) ;
411
412 fclib_problem->info->title = title;
413 fclib_problem->info->description = description;
414 fclib_problem->info->math_info = mathInfo;
415
416
417
418 fclib_problem->spacedim = problem->dimension;
419 fclib_problem->mu = problem->mu;
420 fclib_problem->w = problem->b;
421 fclib_problem->f = problem->q;
422
423 /* only sparse storage */
424 assert(problem->M->matrix2);
425 assert(problem->H->matrix2);
426
427
428 /* only coordinates (triplet) */
429 if(problem->M->matrix2->triplet)
430 {
431 fclib_problem->M = (struct fclib_matrix*)malloc(sizeof(struct fclib_matrix));
432 fclib_problem->M->n = (int) problem->M->matrix2->triplet->n;
433 fclib_problem->M->m = (int) problem->M->matrix2->triplet->m;
434 fclib_problem->M->nzmax= (int) problem->M->matrix2->triplet->nzmax;
435
436 fclib_problem->M->p= (int*) malloc(sizeof(int)*(problem->M->matrix2->triplet->nzmax));
437 csi_to_int(problem->M->matrix2->triplet->p, fclib_problem->M->p,
438 (unsigned) problem->M->matrix2->triplet->nzmax);
439 fclib_problem->M->i= (int*) malloc(sizeof(int)*(problem->M->matrix2->triplet->nzmax));
440 csi_to_int(problem->M->matrix2->triplet->i, fclib_problem->M->i,
441 (unsigned) problem->M->matrix2->triplet->nzmax);
442 fclib_problem->M->x= problem->M->matrix2->triplet->x;
443 fclib_problem->M->nz= (int) problem->M->matrix2->triplet->nz;
444 fclib_problem->M->info=NULL;
445 }
446 else
447 {
448 fprintf(stderr, "globalFrictionContact_fclib_write only implemented for triplet storage.\n");
449 exit(EXIT_FAILURE); ;
450 }
451
452 if(problem->H->matrix2->triplet)
453 {
454 fclib_problem->H = (struct fclib_matrix*)malloc(sizeof(struct fclib_matrix));
455 fclib_problem->H->n = (int) problem->H->matrix2->triplet->n;
456 fclib_problem->H->m = (int) problem->H->matrix2->triplet->m;
457 fclib_problem->H->nzmax= (int) problem->H->matrix2->triplet->nzmax;
458 fclib_problem->H->p= (int*) malloc(sizeof(int)*problem->H->matrix2->triplet->nzmax);
459 csi_to_int(problem->H->matrix2->triplet->p, fclib_problem->H->p,
460 (unsigned) problem->H->matrix2->triplet->nzmax);
461 fclib_problem->H->i= (int*) malloc(sizeof(int)*problem->H->matrix2->triplet->nzmax);
462 csi_to_int(problem->H->matrix2->triplet->i, fclib_problem->H->i,
463 (unsigned) problem->H->matrix2->triplet->nzmax);
464 fclib_problem->H->x= problem->H->matrix2->triplet->x;
465 fclib_problem->H->nz= (int) problem->H->matrix2->triplet->nz;
466 fclib_problem->H->info=NULL;
467 }
468 else
469 {
470 fprintf(stderr, "globalFrictionContact_fclib_write only implemented for triplet storage.\n");
471 exit(EXIT_FAILURE); ;
472 }
473
474 fclib_problem->G = NULL;
475 fclib_problem->b = NULL;
476
477
478
479 DEBUG_PRINT("write in fclib of fclib_problem\n");
480
481 rinfo = fclib_write_global(fclib_problem, path);
482 DEBUG_PRINT("end of write in fclib of fclib_problem\n");
483
484 free(fclib_problem->M->p);
485 free(fclib_problem->M->i);
486 free(fclib_problem->H->p);
487 free(fclib_problem->H->i);
488 free(fclib_problem->H);
489 free(fclib_problem->M);
490 free(fclib_problem->info);
491 free(fclib_problem);
492
493 return rinfo;
494
495 }
496
497
globalRollingFrictionContact_fclib_write(GlobalRollingFrictionContactProblem * problem,char * title,char * description,char * mathInfo,const char * path)498 int globalRollingFrictionContact_fclib_write(
499 GlobalRollingFrictionContactProblem* problem,
500 char * title,
501 char * description,
502 char * mathInfo,
503 const char *path)
504 {
505 int rinfo = 0;
506
507 DEBUG_PRINTF("construction of fclib_problem in %s with title = %s and description = %s\n", path, title, description);
508 if(problem->numberOfContacts == 0)
509 {
510 DEBUG_PRINT("zero contacts");
511 return rinfo;
512
513 }
514 /* globalFrictionContact_display(problem); */
515 /* FILE * file = fopen("toto.dat", "w"); */
516 /* globalFrictionContact_printInFile(problem, file); */
517
518 fclib_global_rolling *fclib_problem;
519 fclib_problem = (fclib_global_rolling*)malloc(sizeof(fclib_global_rolling));
520
521 fclib_problem->info = (struct fclib_info*)malloc(sizeof(struct fclib_info)) ;
522
523 fclib_problem->info->title = title;
524 fclib_problem->info->description = description;
525 fclib_problem->info->math_info = mathInfo;
526
527
528
529 fclib_problem->spacedim = problem->dimension;
530 fclib_problem->mu = problem->mu;
531 fclib_problem->mu_r = problem->mu_r;
532 fclib_problem->w = problem->b;
533 fclib_problem->f = problem->q;
534
535 /* only sparse storage */
536 assert(problem->M->matrix2);
537 assert(problem->H->matrix2);
538
539
540 /* only coordinates (triplet) */
541 if(problem->M->matrix2->triplet)
542 {
543 fclib_problem->M = (struct fclib_matrix*)malloc(sizeof(struct fclib_matrix));
544 fclib_problem->M->n = (int) problem->M->matrix2->triplet->n;
545 fclib_problem->M->m = (int) problem->M->matrix2->triplet->m;
546 fclib_problem->M->nzmax= (int) problem->M->matrix2->triplet->nzmax;
547
548 fclib_problem->M->p= (int*) malloc(sizeof(int)*(problem->M->matrix2->triplet->nzmax));
549 csi_to_int(problem->M->matrix2->triplet->p, fclib_problem->M->p,
550 (unsigned) problem->M->matrix2->triplet->nzmax);
551 fclib_problem->M->i= (int*) malloc(sizeof(int)*(problem->M->matrix2->triplet->nzmax));
552 csi_to_int(problem->M->matrix2->triplet->i, fclib_problem->M->i,
553 (unsigned) problem->M->matrix2->triplet->nzmax);
554 fclib_problem->M->x= problem->M->matrix2->triplet->x;
555 fclib_problem->M->nz= (int) problem->M->matrix2->triplet->nz;
556 fclib_problem->M->info=NULL;
557 }
558 else
559 {
560 fprintf(stderr, "globalRollingFrictionContact_fclib_write only implemented for triplet storage.\n");
561 exit(EXIT_FAILURE); ;
562 }
563
564 if(problem->H->matrix2->triplet)
565 {
566 fclib_problem->H = (struct fclib_matrix*)malloc(sizeof(struct fclib_matrix));
567 fclib_problem->H->n = (int) problem->H->matrix2->triplet->n;
568 fclib_problem->H->m = (int) problem->H->matrix2->triplet->m;
569 fclib_problem->H->nzmax= (int) problem->H->matrix2->triplet->nzmax;
570 fclib_problem->H->p= (int*) malloc(sizeof(int)*problem->H->matrix2->triplet->nzmax);
571 csi_to_int(problem->H->matrix2->triplet->p, fclib_problem->H->p,
572 (unsigned) problem->H->matrix2->triplet->nzmax);
573 fclib_problem->H->i= (int*) malloc(sizeof(int)*problem->H->matrix2->triplet->nzmax);
574 csi_to_int(problem->H->matrix2->triplet->i, fclib_problem->H->i,
575 (unsigned) problem->H->matrix2->triplet->nzmax);
576 fclib_problem->H->x= problem->H->matrix2->triplet->x;
577 fclib_problem->H->nz= (int) problem->H->matrix2->triplet->nz;
578 fclib_problem->H->info=NULL;
579 }
580 else
581 {
582 fprintf(stderr, "globalRollingFrictionContact_fclib_write only implemented for triplet storage.\n");
583 exit(EXIT_FAILURE); ;
584 }
585
586 fclib_problem->G = NULL;
587 fclib_problem->b = NULL;
588
589
590
591 DEBUG_PRINT("write in fclib of fclib_problem\n");
592
593 rinfo = fclib_write_global_rolling(fclib_problem, path);
594 DEBUG_PRINT("end of write in fclib of fclib_problem\n");
595
596 free(fclib_problem->M->p);
597 free(fclib_problem->M->i);
598 free(fclib_problem->H->p);
599 free(fclib_problem->H->i);
600 free(fclib_problem->H);
601 free(fclib_problem->M);
602 free(fclib_problem->info);
603 free(fclib_problem);
604
605 return rinfo;
606
607 }
608
609 #endif
610