Cal_AssembleTerm.cpp

// GetDP - Copyright (C) 1997-2018 P. Dular and C. Geuzaine, University of Liege
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <getdp@onelab.info>.
//
// Contributor(s):
//   Johan Gyselinck
//

#include "ProData.h"
#include "DofData.h"
#include "Message.h"
#include <math.h>

#define SQU(a)     ((a)*(a))

extern struct CurrentData Current ;

static int Warning_Dt = 0, Warning_DtStatic = 0 ;
static int Warning_DtDt = 0, Warning_DtDtStatic = 0, Warning_DtDtFirstOrder = 0 ;

/* ------------------------------------------------------------------------ */
/*  No Time Derivative                                                      */
/* ------------------------------------------------------------------------ */

void Cal_AssembleTerm_NoDt(struct Dof * Equ, struct Dof * Dof, double Val[])
{
  int     k ;
  double  tmp[2] ;

  if(Current.TypeAssembly == ASSEMBLY_SEPARATE){
    if (!Current.DofData->Flag_Init[1]) {
      Current.DofData->Flag_Init[1] = 1 ;
      LinAlg_CreateMatrix(&Current.DofData->M1, &Current.DofData->Solver,
			  Current.DofData->NbrDof, Current.DofData->NbrDof) ;
      LinAlg_CreateVector(&Current.DofData->m1, &Current.DofData->Solver,
			  Current.DofData->NbrDof) ;
      LinAlg_ZeroMatrix(&Current.DofData->M1);
      LinAlg_ZeroVector(&Current.DofData->m1);
      Current.DofData->m1s = List_Create(10, 10, sizeof(gVector));
      for(int i = 0; i < List_Nbr(Current.DofData->TimeFunctionIndex); i++){
        gVector m;
        LinAlg_CreateVector(&m, &Current.DofData->Solver,
                            Current.DofData->NbrDof) ;
        LinAlg_ZeroVector(&m);
        List_Add(Current.DofData->m1s, &m);
      }
    }
    for (k = 0 ; k < Current.NbrHar ; k += 2){
      int incr = (gSCALAR_SIZE == 2) ? k / 2 : k;
      Dof_AssembleInMat(Equ + incr, Dof + incr, Current.NbrHar, &Val[k],
			&Current.DofData->M1, &Current.DofData->m1,
                        Current.DofData->m1s) ;
    }
  }
  else {
    if (Current.NbrHar == 1) {
      switch (Current.TypeTime) {
      case TIME_STATIC :
	Dof_AssembleInMat(Equ, Dof, Current.NbrHar, &Val[0],
			  &Current.DofData->A, &Current.DofData->b) ;
	break ;
      case TIME_THETA :
	tmp[0] = Val[0]*Current.Theta ;
	Dof_AssembleInMat(Equ, Dof, Current.NbrHar, tmp,
			  &Current.DofData->A, &Current.DofData->b) ;
	tmp[0] = Val[0]*(Current.Theta-1.) ;
	Dof_AssembleInVec(Equ, Dof, Current.NbrHar, tmp,
			  Current.DofData->CurrentSolution-1,
			  &(Current.DofData->CurrentSolution-1)->x,