"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "Kernel/Operation_Vector.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.

Operation_Vector.cpp  (getdp-3.4.0-source.tgz)	:	Operation_Vector.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.
#include <map>		#include <map>
#include <string>		#include <string>
#include <vector>		#include <vector>
#include "ProData.h"		#include "ProData.h"
#include "ProParser.h"		#include "ProParser.h"
#include "DofData.h"		#include "DofData.h"
skipping to change at line 42		skipping to change at line 42
*/		*/
// Should preferably be implemented in an object-oriented way		// Should preferably be implemented in an object-oriented way
#define VV_INITIALIZE 0		#define VV_INITIALIZE 0
#define VV_GET_IF_EXISTS 1		#define VV_GET_IF_EXISTS 1
#define VV_GET_AND_ADD_IF_NOT_EXIST 2		#define VV_GET_AND_ADD_IF_NOT_EXIST 2
#define VV_ERASE_IF_EXISTS 3		#define VV_ERASE_IF_EXISTS 3
#define VV_CLEAR_ALL 4		#define VV_CLEAR_ALL 4
int manageVectorMap(int action, char *name, gVector **vector,		int manageVectorMap(int action, char *name, gVector **vector,
struct DofData *DofData_P)		struct DofData *DofData_P)
{		{
// Global map that is never freed		// Global map that is never freed
static std::map<std::string, gVector> vectorMap;		static std::map<std::string, gVector> vectorMap;
if(action == VV_INITIALIZE){ /* Nothing to do*/ }		if(action == VV_INITIALIZE) { /* Nothing to do*/
else if(action == VV_GET_IF_EXISTS)		}
{		else if(action == VV_GET_IF_EXISTS) {
std::map<std::string, gVector>::iterator it = vectorMap.find(name);		std::map<std::string, gVector>::iterator it = vectorMap.find(name);
if(it != vectorMap.end())		if(it != vectorMap.end())
*vector = &it->second;		*vector = &it->second;
else		else
return 1;		return 1;
}		}
else if(action == VV_GET_AND_ADD_IF_NOT_EXIST)		else if(action == VV_GET_AND_ADD_IF_NOT_EXIST) {
{
std::map<std::string, gVector>::iterator it = vectorMap.find(name);		std::map<std::string, gVector>::iterator it = vectorMap.find(name);
if(it != vectorMap.end()){		if(it != vectorMap.end()) { *vector = &it->second; }
*vector = &it->second;		else {
}
else{
gVector n;		gVector n;
LinAlg_CreateVector(&n, &DofData_P->Solver, DofData_P->NbrDof) ;		LinAlg_CreateVector(&n, &DofData_P->Solver, DofData_P->NbrDof);
vectorMap[name] = n;		vectorMap[name] = n;
*vector = &vectorMap[name];		*vector = &vectorMap[name];
}		}
}		}
else if(action == VV_ERASE_IF_EXISTS)		else if(action == VV_ERASE_IF_EXISTS) {
{
std::map<std::string, gVector>::iterator it = vectorMap.find(name);		std::map<std::string, gVector>::iterator it = vectorMap.find(name);
if(it != vectorMap.end())		if(it != vectorMap.end()) {
{
LinAlg_DestroyVector(&it->second);		LinAlg_DestroyVector(&it->second);
vectorMap.erase(name);		vectorMap.erase(name);
}		}
else		else
return 1;		return 1;
}		}
else if(action == VV_CLEAR_ALL)		else if(action == VV_CLEAR_ALL) {
{
for(std::map<std::string, gVector>::iterator it = vectorMap.begin();		for(std::map<std::string, gVector>::iterator it = vectorMap.begin();
it != vectorMap.end(); it++)		it != vectorMap.end(); it++) {
{
LinAlg_DestroyVector(&it->second);		LinAlg_DestroyVector(&it->second);
}		}
vectorMap.clear();		vectorMap.clear();
}		}
else		else {
{
Message::Error("Non valid action for manageVectorMap");		Message::Error("Non valid action for manageVectorMap");
}		}
return 0;		return 0;
}		}
void Operation_CopyVector(struct Operation *Operation_P,		void Operation_CopyVector(struct Operation *Operation_P,
struct DofData *DofData_P)		struct DofData *DofData_P)
{		{
// this global map is never freed for now		// this global map is never freed for now
manageVectorMap(VV_INITIALIZE, NULL, NULL, DofData_P);		manageVectorMap(VV_INITIALIZE, NULL, NULL, DofData_P);
gVector *from = 0, *to = 0, tmp;		gVector *from = 0, *to = 0, tmp;
if(Operation_P->Case.Copy.useList)		if(Operation_P->Case.Copy.useList)
LinAlg_CreateVector(&tmp, &DofData_P->Solver, DofData_P->NbrDof) ;		LinAlg_CreateVector(&tmp, &DofData_P->Solver, DofData_P->NbrDof);
if(Operation_P->Type == OPERATION_COPYSOLUTION){		if(Operation_P->Type == OPERATION_COPYSOLUTION) {
if(DofData_P->CurrentSolution){		if(DofData_P->CurrentSolution) {
if(Operation_P->Case.Copy.from)		if(Operation_P->Case.Copy.from)
to = &DofData_P->CurrentSolution->x;		to = &DofData_P->CurrentSolution->x;
else		else
from = &DofData_P->CurrentSolution->x;		from = &DofData_P->CurrentSolution->x;
}		}
else		else
Message::Error("No current solution available to copy");		Message::Error("No current solution available to copy");
}		}
else if(Operation_P->Type == OPERATION_COPYRHS){		else if(Operation_P->Type == OPERATION_COPYRHS) {
if(Operation_P->Case.Copy.from)		if(Operation_P->Case.Copy.from)
to = &DofData_P->b;		to = &DofData_P->b;
else		else
from = &DofData_P->b;		from = &DofData_P->b;
}		}
else if(Operation_P->Type == OPERATION_COPYRESIDUAL){		else if(Operation_P->Type == OPERATION_COPYRESIDUAL) {
if(Operation_P->Case.Copy.from)		if(Operation_P->Case.Copy.from)
to = &DofData_P->res;		to = &DofData_P->res;
else		else
from = &DofData_P->res;		from = &DofData_P->res;
}		}
else if(Operation_P->Type == OPERATION_COPYINCREMENT){		else if(Operation_P->Type == OPERATION_COPYINCREMENT) {
if(Operation_P->Case.Copy.from)		if(Operation_P->Case.Copy.from)
to = &DofData_P->dx;		to = &DofData_P->dx;
else		else
from = &DofData_P->dx;		from = &DofData_P->dx;
}		}
else{		else {
Message::Error("Copy operation not implemented yet");		Message::Error("Copy operation not implemented yet");
}		}
if(Operation_P->Case.Copy.from){		if(Operation_P->Case.Copy.from) {
if(Operation_P->Case.Copy.useList){		if(Operation_P->Case.Copy.useList) {
Constant *c = Get_ParserConstant(Operation_P->Case.Copy.from);		Constant *c = Get_ParserConstant(Operation_P->Case.Copy.from);
if(c && c->Type == VAR_LISTOFFLOAT){		if(c && c->Type == VAR_LISTOFFLOAT) {
if(List_Nbr(c->Value.List) == DofData_P->NbrDof){		if(List_Nbr(c->Value.List) == DofData_P->NbrDof) {
for(int i = 0; i < List_Nbr(c->Value.List); i++){		for(int i = 0; i < List_Nbr(c->Value.List); i++) {
double d; List_Read(c->Value.List, i, &d);		double d;
		List_Read(c->Value.List, i, &d);
LinAlg_SetDoubleInVector(d, &tmp, i);		LinAlg_SetDoubleInVector(d, &tmp, i);
}		}
LinAlg_AssembleVector(&tmp);		LinAlg_AssembleVector(&tmp);
from = &tmp;		from = &tmp;
}		}
else if(List_Nbr(c->Value.List) == 2 * DofData_P->NbrDof){		else if(List_Nbr(c->Value.List) == 2 * DofData_P->NbrDof) {
for(int i = 0, j = 0; i < List_Nbr(c->Value.List); i += 2, j++){		for(int i = 0, j = 0; i < List_Nbr(c->Value.List); i += 2, j++) {
double d1; List_Read(c->Value.List, i, &d1);		double d1;
double d2; List_Read(c->Value.List, i + 1, &d2);		List_Read(c->Value.List, i, &d1);
		double d2;
		List_Read(c->Value.List, i + 1, &d2);
LinAlg_SetComplexInVector(d1, d2, &tmp, j, j);		LinAlg_SetComplexInVector(d1, d2, &tmp, j, j);
}		}
LinAlg_AssembleVector(&tmp);		LinAlg_AssembleVector(&tmp);
from = &tmp;		from = &tmp;
}		}
else{		else {
Message::Error("Incompatible sizes for vector copy");		Message::Error("Incompatible sizes for vector copy");
}		}
}		}
else if(GetDPNumbers.count(Operation_P->Case.Copy.from)){		else if(GetDPNumbers.count(Operation_P->Case.Copy.from)) {
std::vector<double> &v = GetDPNumbers[Operation_P->Case.Copy.from];		std::vector<double> &v = GetDPNumbers[Operation_P->Case.Copy.from];
if((int)v.size() == DofData_P->NbrDof){		if((int)v.size() == DofData_P->NbrDof) {
for(unsigned int i = 0; i < v.size(); i++){		for(unsigned int i = 0; i < v.size(); i++) {
LinAlg_SetDoubleInVector(v[i], &tmp, i);		LinAlg_SetDoubleInVector(v[i], &tmp, i);
}		}
LinAlg_AssembleVector(&tmp);		LinAlg_AssembleVector(&tmp);
from = &tmp;		from = &tmp;
}		}
else if((int)v.size() == 2 * DofData_P->NbrDof){		else if((int)v.size() == 2 * DofData_P->NbrDof) {
for(unsigned int i = 0; i < v.size(); i += 2){		for(unsigned int i = 0; i < v.size(); i += 2) {
LinAlg_SetComplexInVector(v[i], v[i+1], &tmp, i/2, i/2);		LinAlg_SetComplexInVector(v[i], v[i + 1], &tmp, i / 2, i / 2);
}		}
LinAlg_AssembleVector(&tmp);		LinAlg_AssembleVector(&tmp);
from = &tmp;		from = &tmp;
}		}
else{		else {
Message::Error("Incompatible sizes for vector copy");		Message::Error("Incompatible sizes for vector copy");
}		}
}		}
else{		else {
Message::Error("Non-existant list `%s()' to copy from",		Message::Error("Non-existant list `%s()' to copy from",
Operation_P->Case.Copy.from);		Operation_P->Case.Copy.from);
}		}
}		}
else{		else {
if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.Copy.from, &from, D		if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.Copy.from, &from,
ofData_P))		DofData_P))
Message::Error("Non-existant vector `%s' to copy from",		Message::Error("Non-existant vector `%s' to copy from",
Operation_P->Case.Copy.from);		Operation_P->Case.Copy.from);
}		}
}		}
if(Operation_P->Case.Copy.to){		if(Operation_P->Case.Copy.to) {
if(Operation_P->Case.Copy.useList){		if(Operation_P->Case.Copy.useList) { to = &tmp; }
to = &tmp;		else {
}		manageVectorMap(VV_GET_AND_ADD_IF_NOT_EXIST, Operation_P->Case.Copy.to,
else{		&to, DofData_P);
manageVectorMap(VV_GET_AND_ADD_IF_NOT_EXIST, Operation_P->Case.Copy.to, &t
o, DofData_P);
}		}
}		}
if(from && to){		if(from && to) {
int n1, n2;		int n1, n2;
LinAlg_GetVectorSize(from, &n1);		LinAlg_GetVectorSize(from, &n1);
LinAlg_GetVectorSize(to, &n2);		LinAlg_GetVectorSize(to, &n2);
if(n1 == n2)		if(n1 == n2)
LinAlg_CopyVector(from, to);		LinAlg_CopyVector(from, to);
else		else
Message::Error("Incompatible sizes for vector copy (%d != %d)",		Message::Error("Incompatible sizes for vector copy (%d != %d)", n1, n2);
n1, n2);
}		}
else{		else {
Message::Error("Missing vector for copy");		Message::Error("Missing vector for copy");
}		}
if(Operation_P->Case.Copy.useList){		if(Operation_P->Case.Copy.useList) {
if(Operation_P->Case.Copy.to){		if(Operation_P->Case.Copy.to) {
// create list directly in GetDPNumbers: using parser constants here is		// create list directly in GetDPNumbers: using parser constants here is
// useless since we can never access it		// useless since we can never access it
std::vector<double> v;		std::vector<double> v;
for(int i = 0; i < DofData_P->NbrDof; i++){		for(int i = 0; i < DofData_P->NbrDof; i++) {
gScalar s;		gScalar s;
LinAlg_GetScalarInVector(&s, to, i);		LinAlg_GetScalarInVector(&s, to, i);
if(gSCALAR_SIZE == 2){		if(gSCALAR_SIZE == 2) {
double d1, d2; LinAlg_GetComplexInScalar(&d1, &d2, &s);		double d1, d2;
v.push_back(d1); v.push_back(d2);		LinAlg_GetComplexInScalar(&d1, &d2, &s);
		v.push_back(d1);
		v.push_back(d2);
}		}
else{		else {
double d; LinAlg_GetDoubleInScalar(&d, &s);		double d;
		LinAlg_GetDoubleInScalar(&d, &s);
v.push_back(d);		v.push_back(d);
}		}
}		}
GetDPNumbers[Operation_P->Case.Copy.to] = v;		GetDPNumbers[Operation_P->Case.Copy.to] = v;
}		}
LinAlg_DestroyVector(&tmp) ;		LinAlg_DestroyVector(&tmp);
}		}
}		}
void Operation_AddVector(struct Operation *Operation_P,		void Operation_AddVector(struct Operation *Operation_P,
struct DofData *DofData_P)		struct DofData *DofData_P)
{		{
struct Value Value;		struct Value Value;
Get_ValueOfExpressionByIndex(Operation_P->Case.AddVector.alphaIndex,		Get_ValueOfExpressionByIndex(Operation_P->Case.AddVector.alphaIndex, NULL, 0.,
NULL, 0., 0., 0., &Value);		0., 0., &Value);
double alpha = Value.Val[0];		double alpha = Value.Val[0];
Get_ValueOfExpressionByIndex(Operation_P->Case.AddVector.betaIndex,		Get_ValueOfExpressionByIndex(Operation_P->Case.AddVector.betaIndex, NULL, 0.,
NULL, 0., 0., 0., &Value);		0., 0., &Value);
double beta = Value.Val[0];		double beta = Value.Val[0];
gVector *v1, *v2, *v3;		gVector *v1, *v2, *v3;
// Checking if v1 and v2 exist. If not: error.		// Checking if v1 and v2 exist. If not: error.
if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.AddVector.v1, &v1, DofD		if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.AddVector.v1, &v1,
ata_P))		DofData_P))
Message::Error("Non-existant vector `%s'",		Message::Error("Non-existant vector `%s'", Operation_P->Case.AddVector.v1);
Operation_P->Case.AddVector.v1);		if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.AddVector.v2, &v2,
if(manageVectorMap(VV_GET_IF_EXISTS, Operation_P->Case.AddVector.v2, &v2, DofD		DofData_P))
ata_P))		Message::Error("Non-existant vector `%s'", Operation_P->Case.AddVector.v2);
Message::Error("Non-existant vector `%s'",
Operation_P->Case.AddVector.v2);
// Checking if v3 exists. If not: create it.		// Checking if v3 exists. If not: create it.
manageVectorMap(VV_GET_AND_ADD_IF_NOT_EXIST, Operation_P->Case.AddVector.v3, &		manageVectorMap(VV_GET_AND_ADD_IF_NOT_EXIST, Operation_P->Case.AddVector.v3,
v3, DofData_P);		&v3, DofData_P);
// Check the sizes and perform the operation if OK		// Check the sizes and perform the operation if OK
int n1, n2, n3;		int n1, n2, n3;
LinAlg_GetVectorSize(v1, &n1);		LinAlg_GetVectorSize(v1, &n1);
LinAlg_GetVectorSize(v2, &n2);		LinAlg_GetVectorSize(v2, &n2);
LinAlg_GetVectorSize(v3, &n3);		LinAlg_GetVectorSize(v3, &n3);
if(n1 == n2 && n1 == n3)		if(n1 == n2 && n1 == n3)
LinAlg_AddProdVectorDoubleProdVectorDouble(alpha, v1, beta, v2, v3);		LinAlg_AddProdVectorDoubleProdVectorDouble(alpha, v1, beta, v2, v3);
else		else
Message::Error("Incompatible sizes for vector manipulation (%d, %d and %d)",		Message::Error("Incompatible sizes for vector manipulation (%d, %d and %d)",
n1, n2, n3);		n1, n2, n3);
/* DEBUG (to check the operation is ok)		/* DEBUG (to check the operation is ok)
PetscScalar *a, *b, *c;		PetscScalar *a, *b, *c;
VecGetArray(v1->V, &a); VecGetArray(v2->V, &b); VecGetArray(v3->V, &c);		VecGetArray(v1->V, &a); VecGetArray(v2->V, &b); VecGetArray(v3->V, &c);
for(int i=100; i<110; i++)		for(int i=100; i<110; i++)
{		{
Message::Warning("%g\t%g\t%g", real(a[i]), real(b[i]), real(c[i]));		Message::Warning("%g\t%g\t%g", real(a[i]), real(b[i]), real(c[i]));
}		}
// */		// */
}		}
void Operation_ClearVectors( struct Operation *Operation_P,		void Operation_ClearVectors(struct Operation *Operation_P,
struct DofData *DofData_P)		struct DofData *DofData_P)
{		{
if (List_Nbr(Operation_P->Case.ClearVectors.Names)==0)		if(List_Nbr(Operation_P->Case.ClearVectors.Names) == 0) {
{
Message::Info("ClearVectors: Clear All Vectors");		Message::Info("ClearVectors: Clear All Vectors");
manageVectorMap(VV_CLEAR_ALL, NULL, NULL, DofData_P);		manageVectorMap(VV_CLEAR_ALL, NULL, NULL, DofData_P);
}		}
for(int i = 0; i < List_Nbr(Operation_P->Case.ClearVectors.Names); i++){		for(int i = 0; i < List_Nbr(Operation_P->Case.ClearVectors.Names); i++) {
char *s;		char *s;
List_Read(Operation_P->Case.ClearVectors.Names, i, &s);		List_Read(Operation_P->Case.ClearVectors.Names, i, &s);
if(manageVectorMap(VV_ERASE_IF_EXISTS,s, NULL, DofData_P))		if(manageVectorMap(VV_ERASE_IF_EXISTS, s, NULL, DofData_P))
Message::Info("ClearVectors: Unknown Vector %s", s);		Message::Info("ClearVectors: Unknown Vector %s", s);
else		else
Message::Info("ClearVectors: Clear Vector %s", s);		Message::Info("ClearVectors: Clear Vector %s", s);
}		}
}		}
End of changes. 46 change blocks.
93 lines changed or deleted		85 lines changed or added

---