"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "Kernel/Cal_AssembleTerm.cpp" between
getdp-3.4.0-source.tgz and getdp-3.5.0-source.tgz

About: GetDP is a general finite element solver using mixed elements to discretize de Rham-type complexes in one, two and three dimensions.

[ To the main getdp source changes report ]

Cal_AssembleTerm.cpp  (getdp-3.4.0-source.tgz):Cal_AssembleTerm.cpp  (getdp-3.5.0-source.tgz)
// GetDP - Copyright (C) 1997-2021 P. Dular and C. Geuzaine, University of Liege // GetDP - Copyright (C) 1997-2022 P. Dular and C. Geuzaine, University of Liege
// //
// See the LICENSE.txt file for license information. Please report all // See the LICENSE.txt file for license information. Please report all
// issues on https://gitlab.onelab.info/getdp/getdp/issues. // issues on https://gitlab.onelab.info/getdp/getdp/issues.
// //
// Contributor(s): // Contributor(s):
// Johan Gyselinck // Johan Gyselinck
// //
#include "ProData.h" #include "ProData.h"
#include "DofData.h" #include "DofData.h"
#include "Message.h" #include "Message.h"
#include <math.h> #include <math.h>
#define SQU(a) ((a)*(a)) #define SQU(a) ((a) * (a))
extern struct CurrentData Current ; extern struct CurrentData Current;
static int Warning_Dt = 0, Warning_DtStatic = 0 ; static int Warning_Dt = 0, Warning_DtStatic = 0;
static int Warning_DtDt = 0, Warning_DtDtStatic = 0, Warning_DtDtFirstOrder = 0 static int Warning_DtDt = 0, Warning_DtDtStatic = 0, Warning_DtDtFirstOrder = 0;
;
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* No Time Derivative */ /* No Time Derivative */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_NoDt(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_NoDt(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
double tmp[2] ; double tmp[2];
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[1]) { if(!Current.DofData->Flag_Init[1]) {
Current.DofData->Flag_Init[1] = 1 ; Current.DofData->Flag_Init[1] = 1;
LinAlg_CreateMatrix(&Current.DofData->M1, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M1, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m1, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m1, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M1); LinAlg_ZeroMatrix(&Current.DofData->M1);
LinAlg_ZeroVector(&Current.DofData->m1); LinAlg_ZeroVector(&Current.DofData->m1);
Current.DofData->m1s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m1s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m1s, &m); List_Add(Current.DofData->m1s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2){ for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M1, &Current.DofData->m1, &Current.DofData->M1, &Current.DofData->m1,
Current.DofData->m1s) ; Current.DofData->m1s);
} }
} }
else { else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0], Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0],
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
break ; break;
case TIME_THETA : case TIME_THETA:
tmp[0] = Val[0]*Current.Theta ; tmp[0] = Val[0] * Current.Theta;
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->b);
tmp[0] = Val[0]*(Current.Theta-1.) ; tmp[0] = Val[0] * (Current.Theta - 1.);
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInVec(
Current.DofData->CurrentSolution-1, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
&(Current.DofData->CurrentSolution-1)->x, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&Current.DofData->b) ; break;
break ; case TIME_NEWMARK:
case TIME_NEWMARK : tmp[0] = Val[0] * Current.Beta;
tmp[0] = Val[0]*Current.Beta ; Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->b);
&Current.DofData->A, &Current.DofData->b) ; tmp[0] = Val[0] * (2 * Current.Beta - Current.Gamma - 0.5);
tmp[0] = Val[0]*(2*Current.Beta-Current.Gamma-0.5) ; Dof_AssembleInVec(
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
Current.DofData->CurrentSolution-1, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&(Current.DofData->CurrentSolution-1)->x, tmp[0] = Val[0] * (Current.Gamma - Current.Beta - 0.5);
&Current.DofData->b) ; Dof_AssembleInVec(
tmp[0] = Val[0]*(Current.Gamma-Current.Beta-0.5) ; Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 2,
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, &(Current.DofData->CurrentSolution - 2)->x, &Current.DofData->b);
Current.DofData->CurrentSolution-2, break;
&(Current.DofData->CurrentSolution-2)->x, case TIME_GEAR:
&Current.DofData->b) ; Dof_AssembleInMat(Equ, Dof, Current.NbrHar, Val, &Current.DofData->A,
break ; &Current.DofData->b);
case TIME_GEAR : break;
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, Val,
&Current.DofData->A, &Current.DofData->b) ;
break ;
} }
} }
else { else {
for (k = 0 ; k < Current.NbrHar ; k += 2){ for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
} }
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* First order Time Derivative */ /* First order Time Derivative */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_DtDof(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_DtDof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
double tmp[2] ; double tmp[2];
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[2]) { if(!Current.DofData->Flag_Init[2]) {
Current.DofData->Flag_Init[2] = 1 ; Current.DofData->Flag_Init[2] = 1;
LinAlg_CreateMatrix(&Current.DofData->M2, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M2, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m2, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m2, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M2); LinAlg_ZeroMatrix(&Current.DofData->M2);
LinAlg_ZeroVector(&Current.DofData->m2); LinAlg_ZeroVector(&Current.DofData->m2);
Current.DofData->m2s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m2s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m2s, &m); List_Add(Current.DofData->m2s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2){ for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M2, &Current.DofData->m2, &Current.DofData->M2, &Current.DofData->m2,
Current.DofData->m2s) ; Current.DofData->m2s);
} }
} }
else { else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
if(!Warning_DtStatic){ if(!Warning_DtStatic) {
Message::Info(3, "Discarded DtDof term in static analysis"); Message::Info(3, "Discarded DtDof term in static analysis");
Warning_DtStatic = 1 ; Warning_DtStatic = 1;
} }
break; break;
case TIME_THETA : case TIME_THETA:
tmp[0] = Val[0]/Current.DTime ; tmp[0] = Val[0] / Current.DTime;
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->b);
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInVec(
Current.DofData->CurrentSolution-1, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
&(Current.DofData->CurrentSolution-1)->x, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&Current.DofData->b) ; break;
break ; case TIME_NEWMARK:
case TIME_NEWMARK : tmp[0] = Val[0] * Current.Gamma / Current.DTime;
tmp[0] = Val[0]*Current.Gamma/Current.DTime ; Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->b);
&Current.DofData->A, &Current.DofData->b) ; tmp[0] = Val[0] * (2. * Current.Gamma - 1.) / Current.DTime;
tmp[0] = Val[0]*(2.*Current.Gamma-1.)/Current.DTime ; Dof_AssembleInVec(
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
Current.DofData->CurrentSolution-1, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&(Current.DofData->CurrentSolution-1)->x, tmp[0] = Val[0] * (1. - Current.Gamma) / Current.DTime;
&Current.DofData->b) ; Dof_AssembleInVec(
tmp[0] = Val[0]*(1.-Current.Gamma)/Current.DTime ; Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 2,
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, &(Current.DofData->CurrentSolution - 2)->x, &Current.DofData->b);
Current.DofData->CurrentSolution-2, break;
&(Current.DofData->CurrentSolution-2)->x, case TIME_GEAR:
&Current.DofData->b) ; tmp[0] = Val[0] / (Current.bCorrCoeff * Current.DTime);
break ; Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
case TIME_GEAR : &Current.DofData->b);
tmp[0] = Val[0]/(Current.bCorrCoeff*Current.DTime); for(int i = 0; i < Current.CorrOrder; i++) {
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, tmp[0] = Val[0] * Current.aCorrCoeff[i] /
&Current.DofData->A, &Current.DofData->b) ; (Current.bCorrCoeff * Current.DTime);
for (int i=0; i<Current.CorrOrder; i++) { Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp,
tmp[0] = Val[0]*Current.aCorrCoeff[i]/(Current.bCorrCoeff*Current.DTim Current.DofData->CurrentSolution - 1,
e); &(Current.DofData->CurrentSolution - 1 - i)->x,
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->b);
Current.DofData->CurrentSolution-1,
&(Current.DofData->CurrentSolution-1-i)->x,
&Current.DofData->b) ;
} }
break ; break;
} }
} }
else { else {
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
//printf("====hola idiota=== k %d pul %g \n", k, Current.DofData->Val_Pu // printf("====hola idiota=== k %d pul %g \n", k,
lsation[k/2]); // Current.DofData->Val_Pulsation[k/2]);
tmp[0] = -Val[k+1] * Current.DofData->Val_Pulsation[k/2] ; tmp[0] = -Val[k + 1] * Current.DofData->Val_Pulsation[k / 2];
tmp[1] = Val[k] * Current.DofData->Val_Pulsation[k/2] ; tmp[1] = Val[k] * Current.DofData->Val_Pulsation[k / 2];
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, tmp, Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, tmp,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
} }
} }
} }
} }
void Cal_AssembleTerm_Dt(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_Dt(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
if(!Warning_Dt){ if(!Warning_Dt) {
Message::Warning("Dt not implemented, using DtDof instead"); Message::Warning("Dt not implemented, using DtDof instead");
Warning_Dt = 1 ; Warning_Dt = 1;
} }
Cal_AssembleTerm_DtDof(Equ, Dof, Val); Cal_AssembleTerm_DtDof(Equ, Dof, Val);
} }
/* En preparation ... */ /* En preparation ... */
void Cal_AssembleTerm_DtNL(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_DtNL(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
double tmp[2] ; double tmp[2];
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
Message::Error("DtNL not implemented for separate assembly"); Message::Error("DtNL not implemented for separate assembly");
} }
else { else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
if(!Warning_DtStatic){ if(!Warning_DtStatic) {
Message::Info(3, "Discarded DtDof term in static analysis"); Message::Info(3, "Discarded DtDof term in static analysis");
Warning_DtStatic = 1 ; Warning_DtStatic = 1;
} }
break; break;
case TIME_THETA : case TIME_THETA:
tmp[0] = Val[0]/Current.DTime ; tmp[0] = Val[0] / Current.DTime;
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->b);
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInVec(
Current.DofData->CurrentSolution-1, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
&(Current.DofData->CurrentSolution-1)->x, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&Current.DofData->b) ; break;
break ; case TIME_NEWMARK:
case TIME_NEWMARK : Message::Error(
Message::Error("DtNL not implemented for separate assembly with Newmark s "DtNL not implemented for separate assembly with Newmark scheme");
cheme"); return;
return ; case TIME_GEAR:
case TIME_GEAR :
Message::Error("DtNL not implemented for Gear's method"); Message::Error("DtNL not implemented for Gear's method");
return ; return;
} }
} }
else { else {
Message::Error("DtNL not implemented for separate assembly in harmonic ana Message::Error(
lysis"); "DtNL not implemented for separate assembly in harmonic analysis");
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Second order Time Derivative */ /* Second order Time Derivative */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_DtDtDof(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_DtDtDof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
double tmp[2] ; double tmp[2];
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[3]) { if(!Current.DofData->Flag_Init[3]) {
Current.DofData->Flag_Init[3] = 1 ; Current.DofData->Flag_Init[3] = 1;
LinAlg_CreateMatrix(&Current.DofData->M3, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M3, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m3, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m3, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M3); LinAlg_ZeroMatrix(&Current.DofData->M3);
LinAlg_ZeroVector(&Current.DofData->m3); LinAlg_ZeroVector(&Current.DofData->m3);
Current.DofData->m3s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m3s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m3s, &m); List_Add(Current.DofData->m3s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M3, &Current.DofData->m3, &Current.DofData->M3, &Current.DofData->m3,
Current.DofData->m3s) ; Current.DofData->m3s);
} }
} }
else { else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
if(!Warning_DtDtStatic){ if(!Warning_DtDtStatic) {
Message::Info(3, "Discarded DtDtDof term in static analysis"); Message::Info(3, "Discarded DtDtDof term in static analysis");
Warning_DtDtStatic = 1 ; Warning_DtDtStatic = 1;
} }
break; break;
case TIME_THETA : case TIME_THETA:
if(!Warning_DtDtFirstOrder){ if(!Warning_DtDtFirstOrder) {
Message::Info(3, "Discarded DtDtDof term in Theta time scheme"); Message::Info(3, "Discarded DtDtDof term in Theta time scheme");
Warning_DtDtFirstOrder = 1 ; Warning_DtDtFirstOrder = 1;
} }
break; break;
case TIME_GEAR : case TIME_GEAR:
Message::Error("DtDtDof not implemented for Gear's method"); Message::Error("DtDtDof not implemented for Gear's method");
return ; return;
case TIME_NEWMARK : case TIME_NEWMARK:
tmp[0] = Val[0]/SQU(Current.DTime) ; tmp[0] = Val[0] / SQU(Current.DTime);
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->A,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->b);
tmp[0] = 2*Val[0]/SQU(Current.DTime) ; tmp[0] = 2 * Val[0] / SQU(Current.DTime);
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInVec(
Current.DofData->CurrentSolution-1, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 1,
&(Current.DofData->CurrentSolution-1)->x, &(Current.DofData->CurrentSolution - 1)->x, &Current.DofData->b);
&Current.DofData->b) ; tmp[0] = -Val[0] / SQU(Current.DTime);
tmp[0] = -Val[0]/SQU(Current.DTime) ; Dof_AssembleInVec(
Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp, Equ, Dof, Current.NbrHar, tmp, Current.DofData->CurrentSolution - 2,
Current.DofData->CurrentSolution-2, &(Current.DofData->CurrentSolution - 2)->x, &Current.DofData->b);
&(Current.DofData->CurrentSolution-2)->x, break;
&Current.DofData->b) ;
break ;
} }
} }
else { else {
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
tmp[0] = - Val[k] * SQU(Current.DofData->Val_Pulsation[k/2]) ; tmp[0] = -Val[k] * SQU(Current.DofData->Val_Pulsation[k / 2]);
tmp[1] = - Val[k+1] * SQU(Current.DofData->Val_Pulsation[k/2]) ; tmp[1] = -Val[k + 1] * SQU(Current.DofData->Val_Pulsation[k / 2]);
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, tmp, Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, tmp,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
} }
} }
} }
} }
void Cal_AssembleTerm_DtDt(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_DtDt(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
if(!Warning_DtDt){ if(!Warning_DtDt) {
Message::Warning("DtDt not implemented, using DtDtDof instead"); Message::Warning("DtDt not implemented, using DtDtDof instead");
Warning_DtDt = 1 ; Warning_DtDt = 1;
} }
Cal_AssembleTerm_DtDtDof(Equ, Dof, Val); Cal_AssembleTerm_DtDtDof(Equ, Dof, Val);
} }
/* ------------------------------------------------- */ /* ------------------------------------------------- */
/* higher order Time Derivative for Polynomial EVP */ /* higher order Time Derivative for Polynomial EVP */
/* ------------------------------------------------- */ /* ------------------------------------------------- */
void Cal_AssembleTerm_DtDtDtDof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_DtDtDtDof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[4]) { if(!Current.DofData->Flag_Init[4]) {
Current.DofData->Flag_Init[4] = 1 ; Current.DofData->Flag_Init[4] = 1;
LinAlg_CreateMatrix(&Current.DofData->M4, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M4, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m4, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m4, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M4); LinAlg_ZeroMatrix(&Current.DofData->M4);
LinAlg_ZeroVector(&Current.DofData->m4); LinAlg_ZeroVector(&Current.DofData->m4);
Current.DofData->m4s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m4s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m4s, &m); List_Add(Current.DofData->m4s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M4, &Current.DofData->m4, &Current.DofData->M4, &Current.DofData->m4,
Current.DofData->m4s) ; Current.DofData->m4s);
} }
} }
else { else {
Message::Error("DtDtDtDof only available with GenerateSeparate"); Message::Error("DtDtDtDof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_DtDtDtDtDof(struct Dof * Equ, struct Dof * Dof, double Val void Cal_AssembleTerm_DtDtDtDtDof(struct Dof *Equ, struct Dof *Dof,
[]) double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[5]) { if(!Current.DofData->Flag_Init[5]) {
Current.DofData->Flag_Init[5] = 1 ; Current.DofData->Flag_Init[5] = 1;
LinAlg_CreateMatrix(&Current.DofData->M5, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M5, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m5, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m5, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M5); LinAlg_ZeroMatrix(&Current.DofData->M5);
LinAlg_ZeroVector(&Current.DofData->m5); LinAlg_ZeroVector(&Current.DofData->m5);
Current.DofData->m5s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m5s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m5s, &m); List_Add(Current.DofData->m5s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M5, &Current.DofData->m5, &Current.DofData->M5, &Current.DofData->m5,
Current.DofData->m5s) ; Current.DofData->m5s);
} }
} }
else { else {
Message::Error("DtDtDtDtDof only available with GenerateSeparate"); Message::Error("DtDtDtDtDof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_DtDtDtDtDtDof(struct Dof * Equ, struct Dof * Dof, double V void Cal_AssembleTerm_DtDtDtDtDtDof(struct Dof *Equ, struct Dof *Dof,
al[]) double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[6]) { if(!Current.DofData->Flag_Init[6]) {
Current.DofData->Flag_Init[6] = 1 ; Current.DofData->Flag_Init[6] = 1;
LinAlg_CreateMatrix(&Current.DofData->M6, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M6, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m6, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m6, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M6); LinAlg_ZeroMatrix(&Current.DofData->M6);
LinAlg_ZeroVector(&Current.DofData->m6); LinAlg_ZeroVector(&Current.DofData->m6);
Current.DofData->m6s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m6s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m6s, &m); List_Add(Current.DofData->m6s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M6, &Current.DofData->m6, &Current.DofData->M6, &Current.DofData->m6,
Current.DofData->m6s) ; Current.DofData->m6s);
} }
} }
else { else {
Message::Error("DtDtDtDtDtDof only available with GenerateSeparate"); Message::Error("DtDtDtDtDtDof only available with GenerateSeparate");
return ; return;
} }
} }
// nleigchange // nleigchange
void Cal_AssembleTerm_NLEig1Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig1Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[7]) { if(!Current.DofData->Flag_Init[7]) {
Current.DofData->Flag_Init[7] = 1 ; Current.DofData->Flag_Init[7] = 1;
LinAlg_CreateMatrix(&Current.DofData->M7, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M7, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m7, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m7, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M7); LinAlg_ZeroMatrix(&Current.DofData->M7);
LinAlg_ZeroVector(&Current.DofData->m7); LinAlg_ZeroVector(&Current.DofData->m7);
Current.DofData->m7s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m7s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m7s, &m); List_Add(Current.DofData->m7s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M7, &Current.DofData->m7, &Current.DofData->M7, &Current.DofData->m7,
Current.DofData->m7s) ; Current.DofData->m7s);
} }
} }
else { else {
Message::Error("NLEig1Dof only available with GenerateSeparate"); Message::Error("NLEig1Dof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_NLEig2Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig2Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[6]) { if(!Current.DofData->Flag_Init[6]) {
Current.DofData->Flag_Init[6] = 1 ; Current.DofData->Flag_Init[6] = 1;
LinAlg_CreateMatrix(&Current.DofData->M6, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M6, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m6, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m6, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M6); LinAlg_ZeroMatrix(&Current.DofData->M6);
LinAlg_ZeroVector(&Current.DofData->m6); LinAlg_ZeroVector(&Current.DofData->m6);
Current.DofData->m6s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m6s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m6s, &m); List_Add(Current.DofData->m6s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M6, &Current.DofData->m6, &Current.DofData->M6, &Current.DofData->m6,
Current.DofData->m6s) ; Current.DofData->m6s);
} }
} }
else { else {
Message::Error("NLEig3Dof only available with GenerateSeparate"); Message::Error("NLEig3Dof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_NLEig3Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig3Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[5]) { if(!Current.DofData->Flag_Init[5]) {
Current.DofData->Flag_Init[5] = 1 ; Current.DofData->Flag_Init[5] = 1;
LinAlg_CreateMatrix(&Current.DofData->M5, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M5, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m5, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m5, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M5); LinAlg_ZeroMatrix(&Current.DofData->M5);
LinAlg_ZeroVector(&Current.DofData->m5); LinAlg_ZeroVector(&Current.DofData->m5);
Current.DofData->m5s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m5s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m5s, &m); List_Add(Current.DofData->m5s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M5, &Current.DofData->m5, &Current.DofData->M5, &Current.DofData->m5,
Current.DofData->m5s) ; Current.DofData->m5s);
} }
} }
else { else {
Message::Error("NLEig3Dof only available with GenerateSeparate"); Message::Error("NLEig3Dof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_NLEig4Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig4Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[4]) { if(!Current.DofData->Flag_Init[4]) {
Current.DofData->Flag_Init[4] = 1 ; Current.DofData->Flag_Init[4] = 1;
LinAlg_CreateMatrix(&Current.DofData->M4, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M4, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m4, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m4, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M4); LinAlg_ZeroMatrix(&Current.DofData->M4);
LinAlg_ZeroVector(&Current.DofData->m4); LinAlg_ZeroVector(&Current.DofData->m4);
Current.DofData->m4s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m4s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m4s, &m); List_Add(Current.DofData->m4s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M4, &Current.DofData->m4, &Current.DofData->M4, &Current.DofData->m4,
Current.DofData->m4s) ; Current.DofData->m4s);
} }
} }
else { else {
Message::Error("NLEig4Dof only available with GenerateSeparate"); Message::Error("NLEig4Dof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_NLEig5Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig5Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[3]) { if(!Current.DofData->Flag_Init[3]) {
Current.DofData->Flag_Init[3] = 1 ; Current.DofData->Flag_Init[3] = 1;
LinAlg_CreateMatrix(&Current.DofData->M3, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M3, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m3, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m3, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M3); LinAlg_ZeroMatrix(&Current.DofData->M3);
LinAlg_ZeroVector(&Current.DofData->m3); LinAlg_ZeroVector(&Current.DofData->m3);
Current.DofData->m3s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m3s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m3s, &m); List_Add(Current.DofData->m3s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M3, &Current.DofData->m3, &Current.DofData->M3, &Current.DofData->m3,
Current.DofData->m3s) ; Current.DofData->m3s);
} }
} }
else { else {
Message::Error("NLEig5Dof only available with GenerateSeparate"); Message::Error("NLEig5Dof only available with GenerateSeparate");
return ; return;
} }
} }
void Cal_AssembleTerm_NLEig6Dof(struct Dof * Equ, struct Dof * Dof, double Val[] ) void Cal_AssembleTerm_NLEig6Dof(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[2]) { if(!Current.DofData->Flag_Init[2]) {
Current.DofData->Flag_Init[2] = 1 ; Current.DofData->Flag_Init[2] = 1;
LinAlg_CreateMatrix(&Current.DofData->M2, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M2, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m2, &Current.DofDat LinAlg_CreateVector(&Current.DofData->m2, &Current.DofData->Solver,
a->Solver, Current.DofData->NbrDof);
Current.DofData->NbrDof) ;
LinAlg_ZeroMatrix(&Current.DofData->M2); LinAlg_ZeroMatrix(&Current.DofData->M2);
LinAlg_ZeroVector(&Current.DofData->m2); LinAlg_ZeroVector(&Current.DofData->m2);
Current.DofData->m2s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m2s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m2s, &m); List_Add(Current.DofData->m2s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M2, &Current.DofData->m2, &Current.DofData->M2, &Current.DofData->m2,
Current.DofData->m2s) ; Current.DofData->m2s);
} }
} }
else { else {
Message::Error("NLEig6Dof only available with GenerateSeparate"); Message::Error("NLEig6Dof only available with GenerateSeparate");
return ; return;
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Jacobian NonLinear */ /* Jacobian NonLinear */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_JacNL(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_JacNL(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
Message::Error("JacNL not implemented for separate assembly"); Message::Error("JacNL not implemented for separate assembly");
} }
else{ else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
case TIME_THETA : case TIME_THETA:
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0],
&Current.DofData->Jac, NULL);
break;
case TIME_GEAR:
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0], Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0],
&Current.DofData->Jac, NULL) ; &Current.DofData->Jac, NULL);
break ; break;
case TIME_GEAR : case TIME_NEWMARK:
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0], Message::Error("JacNL not implemented for Newmark's method");
&Current.DofData->Jac, NULL) ; return;
break ;
case TIME_NEWMARK :
Message::Error("JacNL not implemented for Newmark's method");
return ;
} }
} }
else { else {
for (k = 0 ; k < Current.NbrHar ; k += 2){ for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->Jac, NULL) ; &Current.DofData->Jac, NULL);
} }
} }
} }
} }
void Cal_AssembleTerm_DtDofJacNL(struct Dof * Equ, struct Dof * Dof, double Val[ ]) void Cal_AssembleTerm_DtDofJacNL(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
double tmp[2] ; double tmp[2];
if(Current.TypeAssembly == ASSEMBLY_SEPARATE) if(Current.TypeAssembly == ASSEMBLY_SEPARATE)
Message::Error("DtDofJacNL not implemented for separate assembly"); Message::Error("DtDofJacNL not implemented for separate assembly");
else { else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
if(!Warning_DtStatic){ if(!Warning_DtStatic) {
Message::Info(3, "Discarded DtDofJacNL term in static analysis"); Message::Info(3, "Discarded DtDofJacNL term in static analysis");
Warning_DtStatic = 1 ; Warning_DtStatic = 1;
} }
break; break;
case TIME_THETA : case TIME_THETA:
if ( fabs(Current.Theta - 1.0) > 1e-3 ){ if(fabs(Current.Theta - 1.0) > 1e-3) {
Message::Error("Theta method not implemented for nonlinear problems wh Message::Error(
en " "Theta method not implemented for nonlinear problems when "
"Theta != 1.0"); "Theta != 1.0");
return; return;
} }
tmp[0] = Val[0]/Current.DTime; tmp[0] = Val[0] / Current.DTime;
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->Jac,
&Current.DofData->Jac, NULL); NULL);
break ; break;
case TIME_NEWMARK : case TIME_NEWMARK:
Message::Error("DtDofJacNL not implemented for Newmark scheme"); Message::Error("DtDofJacNL not implemented for Newmark scheme");
return ; return;
case TIME_GEAR : case TIME_GEAR:
tmp[0] = Val[0] / (Current.bCorrCoeff * Current.DTime); tmp[0] = Val[0] / (Current.bCorrCoeff * Current.DTime);
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp, &Current.DofData->Jac,
&Current.DofData->Jac, NULL); NULL);
break ; break;
} }
} }
else else
Message::Error("DtDofJacNL not implemented for harmonic analysis"); Message::Error("DtDofJacNL not implemented for harmonic analysis");
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Never Time Derivative (provisoire mais tres important ... Patrick) */ /* Never Time Derivative (provisoire mais tres important ... Patrick) */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_NeverDt(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_NeverDt(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
int k ; int k;
if(Current.TypeAssembly == ASSEMBLY_SEPARATE){ if(Current.TypeAssembly == ASSEMBLY_SEPARATE) {
if (!Current.DofData->Flag_Init[1]) { if(!Current.DofData->Flag_Init[1]) {
Current.DofData->Flag_Init[1] = 1 ; Current.DofData->Flag_Init[1] = 1;
LinAlg_CreateMatrix(&Current.DofData->M1, &Current.DofData->Solver, LinAlg_CreateMatrix(&Current.DofData->M1, &Current.DofData->Solver,
Current.DofData->NbrDof, Current.DofData->NbrDof) ; Current.DofData->NbrDof, Current.DofData->NbrDof);
LinAlg_CreateVector(&Current.DofData->m1, &Current.DofData->Solver, LinAlg_CreateVector(&Current.DofData->m1, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroMatrix(&Current.DofData->M1); LinAlg_ZeroMatrix(&Current.DofData->M1);
LinAlg_ZeroVector(&Current.DofData->m1); LinAlg_ZeroVector(&Current.DofData->m1);
Current.DofData->m1s = List_Create(10, 10, sizeof(gVector)); Current.DofData->m1s = List_Create(10, 10, sizeof(gVector));
for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){ for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++) {
gVector m; gVector m;
LinAlg_CreateVector(&m, &Current.DofData->Solver, LinAlg_CreateVector(&m, &Current.DofData->Solver,
Current.DofData->NbrDof) ; Current.DofData->NbrDof);
LinAlg_ZeroVector(&m); LinAlg_ZeroVector(&m);
List_Add(Current.DofData->m1s, &m); List_Add(Current.DofData->m1s, &m);
} }
} }
for (k = 0 ; k < Current.NbrHar ; k += 2){ for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->M1, &Current.DofData->m1, &Current.DofData->M1, &Current.DofData->m1,
Current.DofData->m1s) ; Current.DofData->m1s);
} }
} }
else{ else {
if (Current.NbrHar == 1) { if(Current.NbrHar == 1) {
switch (Current.TypeTime) { switch(Current.TypeTime) {
case TIME_STATIC : case TIME_STATIC:
case TIME_THETA : case TIME_THETA:
case TIME_NEWMARK : case TIME_NEWMARK:
case TIME_GEAR : case TIME_GEAR:
Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0], Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0],
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
break ; break;
} }
} }
else { else {
for (k = 0 ; k < Current.NbrHar ; k += 2) { for(k = 0; k < Current.NbrHar; k += 2) {
int incr = (gSCALAR_SIZE == 2) ? k / 2 : k; int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k], Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->A, &Current.DofData->b);
} }
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Multi-Harmonic with movement */ /* Multi-Harmonic with movement */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Cal_AssembleTerm_MHMoving(struct Dof * Equ, struct Dof * Dof, double Val[]) void Cal_AssembleTerm_MHMoving(struct Dof *Equ, struct Dof *Dof, double Val[])
{ {
// MHMoving_assemblyType = 1 => Use current system A,b // MHMoving_assemblyType = 1 => Use current system A,b
// MHMoving_assemblyType = 2 => Use dedicated system A_MH_Moving, b_MH_Moving // MHMoving_assemblyType = 2 => Use dedicated system A_MH_Moving, b_MH_Moving
// MHMoving_assemblyType = 3 => Look for unknowns and constraints in Moving Gr // MHMoving_assemblyType = 3 => Look for unknowns and constraints in Moving
oup // Group
extern int MHMoving_assemblyType ; extern int MHMoving_assemblyType;
extern double ** MH_Moving_Matrix ; extern double **MH_Moving_Matrix;
extern Tree_T * DofTree_MH_moving ; extern Tree_T *DofTree_MH_moving;
// FIXME: this cannot work in complex arithmetic: AssembleInMat will // FIXME: this cannot work in complex arithmetic: AssembleInMat will
// need to assemble both real and imaginary parts at once -- See // need to assemble both real and imaginary parts at once -- See
// Cal_AssembleTerm for an example // Cal_AssembleTerm for an example
if(MHMoving_assemblyType==1){ if(MHMoving_assemblyType == 1) {
for (int k = 0 ; k < Current.NbrHar ; k++) for(int k = 0; k < Current.NbrHar; k++)
for (int l = 0 ; l < Current.NbrHar ; l++) { for(int l = 0; l < Current.NbrHar; l++) {
double tmp = Val[0] * MH_Moving_Matrix[k][l] ; double tmp = Val[0] * MH_Moving_Matrix[k][l];
// if (k==l) // if (k==l)
Dof_AssembleInMat(Equ+k, Dof+l, 1, &tmp, Dof_AssembleInMat(Equ + k, Dof + l, 1, &tmp, &Current.DofData->A,
&Current.DofData->A, &Current.DofData->b) ; &Current.DofData->b);
} }
} }
if(MHMoving_assemblyType==2){ if(MHMoving_assemblyType == 2) {
for (int k = 0 ; k < Current.NbrHar ; k++) for(int k = 0; k < Current.NbrHar; k++)
for (int l = 0 ; l < Current.NbrHar ; l++) { for(int l = 0; l < Current.NbrHar; l++) {
double tmp = Val[0] * MH_Moving_Matrix[k][l] ; double tmp = Val[0] * MH_Moving_Matrix[k][l];
// if (k==l) // if (k==l)
Dof_AssembleInMat(Equ+k, Dof+l, 1, &tmp, &Current.DofData->A_MH_moving, Dof_AssembleInMat(Equ + k, Dof + l, 1, &tmp,
NULL); &Current.DofData->A_MH_moving, NULL);
// &Current.DofData->A_MH_moving, &Current.DofData->b_ // &Current.DofData->A_MH_moving, &Current.DofData->b_MH_moving) ;
MH_moving) ;
} }
} }
if(MHMoving_assemblyType==3){ if(MHMoving_assemblyType == 3) {
if (Dof->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Dof)) if(Dof->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Dof))
Tree_Add(DofTree_MH_moving,Dof) ; Tree_Add(DofTree_MH_moving, Dof);
else if (Dof->Type == DOF_LINK && !Tree_PQuery(DofTree_MH_moving, Dof->Case. else if(Dof->Type == DOF_LINK &&
Link.Dof)) !Tree_PQuery(DofTree_MH_moving, Dof->Case.Link.Dof))
Tree_Add(DofTree_MH_moving,Dof->Case.Link.Dof) ; Tree_Add(DofTree_MH_moving, Dof->Case.Link.Dof);
if (Equ->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Equ)) if(Equ->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Equ))
Tree_Add(DofTree_MH_moving,Equ) ; Tree_Add(DofTree_MH_moving, Equ);
else if (Equ->Type == DOF_LINK && !Tree_PQuery(DofTree_MH_moving, Equ->Case. else if(Equ->Type == DOF_LINK &&
Link.Dof)) !Tree_PQuery(DofTree_MH_moving, Equ->Case.Link.Dof))
Tree_Add(DofTree_MH_moving,Equ->Case.Link.Dof) ; Tree_Add(DofTree_MH_moving, Equ->Case.Link.Dof);
} }
} }
/*
void Cal_AssembleTerm_MH_Moving_simple(struct Dof * Equ, struct Dof * Dof, doubl
e Val[])
{
int k, l ;
double tmp ;
extern double ** MH_Moving_Matrix ;
for (k = 0 ; k < Current.NbrHar ; k++)
for (l = 0 ; l < Current.NbrHar ; l++) {
tmp = Val[0] * MH_Moving_Matrix[k][l] ;
// if (k==l)
Dof_AssembleInMat(Equ+k, Dof+l, 1, &tmp,
&Current.DofData->A, &Current.DofData->b) ;
}
}
void Cal_AssembleTerm_MH_Moving_separate(struct Dof * Equ, struct Dof * Dof, dou
ble Val[])
{
int k, l ;
double tmp ;
extern double ** MH_Moving_Matrix ;
for (k = 0 ; k < Current.NbrHar ; k++)
for (l = 0 ; l < Current.NbrHar ; l++) {
tmp = Val[0] * MH_Moving_Matrix[k][l] ;
// if (k==l)
Dof_AssembleInMat(Equ+k, Dof+l, 1, &tmp,
&Current.DofData->A_MH_moving, &Current.DofData->b_MH_mov
ing) ;
}
}
void Cal_AssembleTerm_MH_Moving_probe(struct Dof * Equ, struct Dof * Dof, double
Val[])
{
extern Tree_T * DofTree_MH_moving ;
if (Dof->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Dof))
Tree_Add(DofTree_MH_moving,Dof) ;
else if (Dof->Type == DOF_LINK && !Tree_PQuery(DofTree_MH_moving, Dof->Case.Li
nk.Dof))
Tree_Add(DofTree_MH_moving,Dof->Case.Link.Dof) ;
if (Equ->Type == DOF_UNKNOWN && !Tree_PQuery(DofTree_MH_moving, Equ))
Tree_Add(DofTree_MH_moving,Equ) ;
else if (Equ->Type == DOF_LINK && !Tree_PQuery(DofTree_MH_moving, Equ->Case.Li
nk.Dof))
Tree_Add(DofTree_MH_moving,Equ->Case.Link.Dof) ;
}
*/
 End of changes. 163 change blocks. 
436 lines changed or deleted 418 lines changed or added
To the C++ Programs section of the getdp source changes report or the top of this page
Mainly generated by pkgdiff using rfcdiff
Home  |  About  |  Features  |  All  |  Newest  |  Dox  |  Diffs  |  RSS Feeds  |  Screenshots  |  Comments  |  Imprint  |  Privacy  |  HTTP(S)