getdp issueshttps://gitlab.onelab.info/getdp/getdp/issues2019-07-09T17:35:40Zhttps://gitlab.onelab.info/getdp/getdp/issues/57Link Constraint for Non-Identical Regions2019-07-09T17:35:40ZBruno de Sousa AlvesLink Constraint for Non-Identical RegionsDear all,
I hope this is the correct place to post a question like that.
I’ve been working in a combined 1D and 2D finite element (FE) model in GetDP with a magnetic field (H) formulation and I’d like to connect the intensity of magnetic field over an edge in 2D with the boundary of a 1D FE system (see figure below). For this purpose, I’ve created two independent “FunctionSpaces” and I tried to apply a constraint of type “Link” to connect the related Degrees of Freedom. However, as each edge in the 2D system is defined by two nodes, the elements in the 2D and 1D regions are not geometrically identical (edges in 2D and nodes in 1D) and I couldn’t apply this type of constraint directly.
![Link2](/uploads/8bad1305cb113bf73c78c7abea7d8360/Link2.jpg)
Thus, I was wondering if there is any way to connect them using a constraint “Link”. Maybe with the optional “FunctionRef” function, but I couldn’t figure out how it works. If possible, could you please give me an example of the application of the “Link” with “FunctionRef”?
Thank you and best regards,
Bruno AlvesDear all,
I hope this is the correct place to post a question like that.
I’ve been working in a combined 1D and 2D finite element (FE) model in GetDP with a magnetic field (H) formulation and I’d like to connect the intensity of magnetic field over an edge in 2D with the boundary of a 1D FE system (see figure below). For this purpose, I’ve created two independent “FunctionSpaces” and I tried to apply a constraint of type “Link” to connect the related Degrees of Freedom. However, as each edge in the 2D system is defined by two nodes, the elements in the 2D and 1D regions are not geometrically identical (edges in 2D and nodes in 1D) and I couldn’t apply this type of constraint directly.
![Link2](/uploads/8bad1305cb113bf73c78c7abea7d8360/Link2.jpg)
Thus, I was wondering if there is any way to connect them using a constraint “Link”. Maybe with the optional “FunctionRef” function, but I couldn’t figure out how it works. If possible, could you please give me an example of the application of the “Link” with “FunctionRef”?
Thank you and best regards,
Bruno Alveshttps://gitlab.onelab.info/getdp/getdp/issues/55Computing rotor speed and position2019-05-20T15:17:23ZCássio Krugerkrugercassio@gmail.comComputing rotor speed and positionHello guys! I've been work on a magnetic gear project (https://github.com/CassioKruger/PDD-pure) and I'm having trouble to compute the velocity of the moving bands. I took the "machine_magstadyn_a.pro" file to base my project and added some modifications to use it with 2 rotors instead of 1.
There is a gear ratio between both rotors, so I made this in my .pro file:
```
// pdd gear ratio:
// pH*wH + pL*wL = nP*wP
// wH is the speed of the inner rotor
// wL is the speed of the outer rotor
// wP is the speed of the modulators
// When one of the three parts of the gear is stationary, there will be a constant relation or gear ratio
// between the speeds of other two parts.
// considering that the outer rotor is stationary, the gear ratio becomes:
// -> pH*wH = nP*wP
// -> Gr = pH/nP = wP/wH
// -> gear ratio = nbr of poles at rotor 1 / nbr of modulators
// in this case, the nbr of modulators is equal to the nbr of poles at the outer rotor, so:
gear_ratio = NbrPolesInModel/NbrSectStatorMag;
delta_theta[] = delta_theta_deg * deg2rad ; //rotor 1 step
delta_theta2[] = delta_theta[] * gear_ratio ; //rotor 2 step
```
those "delta_theta[]" should be the step of each rotor to use with "ChangeOfCoordinates" inside the timeloop of the "machine_magstadyn_a_2rotors.pro" file, like this:
```
ChangeOfCoordinates[ NodesOf[Rotor_Moving], RotatePZ[delta_theta[]]];
ChangeOfCoordinates[ NodesOf[Rotor2_Moving], RotatePZ[delta_theta2[]]];
```
and that's working nice, the gear ration is clearly seen at the simulation. But the thing that I'm couldn't figure out is how to compute the correct speed and position of each rotor, so I can show it with graphs and prove that my model is working. I've made this:
-First, I've declared a "DomainKin2", because now there is 2 parts moving
-Then, inside the "Constraint{}" I made this:
```
Constraint{
...
...
...
//Kinetics
{ Name CurrentPosition ; // [m]
Case {
{ Region DomainKin ; Type Init ; Value ($PreviousPosition = 0) ; }
}
}
{ Name CurrentVelocity ; // [rad/s]
Case {
{ Region DomainKin ; Type Init ; Value ($PreviousVelocity = 0) ; }
}
}
//Kinetics - MOVING BAND 2
{ Name CurrentPosition2 ; // [m]
Case {
{ Region DomainKin2 ; Type Init ; Value ($PreviousPosition2 = 0) ; }
}
}
{ Name CurrentVelocity2 ; // [rad/s]
Case {
{ Region DomainKin2 ; Type Init ; Value ($PreviousVelocity2 = 0) ; }
}
}
}
```
-After that, I did the same to "FunctionSpace{}":
```
FunctionSpace{
...
...
...
// For mechanical equation
{ Name Position ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef pr ; Function BF_Region ;
Support DomainKin ; Entity DomainKin ; }
}
GlobalQuantity {
{ Name P ; Type AliasOf ; NameOfCoef pr ; }
}
Constraint {
{ NameOfCoef P ; EntityType Region ; NameOfConstraint CurrentPosition ; }
}
}
{ Name Velocity ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef vr ; Function BF_Region ;
Support DomainKin ; Entity DomainKin ; } }
GlobalQuantity {
{ Name V ; Type AliasOf ; NameOfCoef vr ; }
}
Constraint {
{ NameOfCoef V ; EntityType Region ; NameOfConstraint CurrentVelocity ; }
}
}
// For mechanical equation - MOVING BAND 2
{ Name Position2 ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef pr2 ; Function BF_Region ;
Support DomainKin2 ; Entity DomainKin2 ; }
}
GlobalQuantity {
{ Name P2 ; Type AliasOf ; NameOfCoef pr2 ; }
}
Constraint {
{ NameOfCoef P2 ; EntityType Region ; NameOfConstraint CurrentPosition2 ; }
}
}
{ Name Velocity2 ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef vr2 ; Function BF_Region ;
Support DomainKin2 ; Entity DomainKin2 ; } }
GlobalQuantity {
{ Name V2 ; Type AliasOf ; NameOfCoef vr2 ; }
}
Constraint {
{ NameOfCoef V2 ; EntityType Region ; NameOfConstraint CurrentVelocity2 ; }
}
}
}
```
-In the "Formulation{}", I did this:
```
Formulation{
...
...
...
// Mechanics
{ Name Mechanical ; Type FemEquation ;
Quantity {
{ Name V ; Type Global ; NameOfSpace Velocity [V] ; } // velocity
{ Name P ; Type Global ; NameOfSpace Position [P] ; } // position
{ Name V2 ; Type Global ; NameOfSpace Velocity2 [V2] ; } // velocity MB2
{ Name P2 ; Type Global ; NameOfSpace Position2 [P2] ; } // position MB2
}
Equation {
GlobalTerm { DtDof [ /*Inertia **/ Dof{V} , {V} ] ; In DomainKin ; }
//GlobalTerm { [ Friction[] * Dof{V} , {V} ] ; In DomainKin ; }
GlobalTerm { [ Torque_mec[], {V} ] ; In DomainKin ; }
GlobalTerm { [ Torque_mag[] , {V} ] ; In DomainKin ; }
GlobalTerm { DtDof [ Dof{P} , {P} ] ; In DomainKin ; }
GlobalTerm { [-Dof{V} , {P} ] ; In DomainKin ; }
//---------------------------------------------------------------//
GlobalTerm { DtDof [ /*Inertia*0.7 **/ Dof{V2} , {V2} ] ; In DomainKin2 ; }
//GlobalTerm { [ Friction[] * Dof{V2} , {V2} ] ; In DomainKin2 ; }
GlobalTerm { [ Torque_mec[] , {V2} ] ; In DomainKin2 ; }
GlobalTerm { [ Torque_mag[] , {V2} ] ; In DomainKin2 ; }
GlobalTerm { DtDof [ Dof{P2} , {P2} ] ; In DomainKin2 ; }
GlobalTerm { [-Dof{V2} , {P2} ] ; In DomainKin2 ; }
}
}
}
```
-And, finally, in the "PostProcessing{}" I did this:
```
PostProcessing{
...
...
...
{ Name Mechanical ; NameOfFormulation Mechanical ;
PostQuantity {
{ Name P ; Value { Term { [ {P} ] ; In DomainKin ; } } } // Position [rad]
{ Name Pdeg ; Value { Term { [ {P}*180/Pi ] ; In DomainKin ; } } } // Position [deg]
{ Name V ; Value { Term { [ {V} ] ; In DomainKin ; } } } // Velocity [rad/s]
{ Name Vrpm ; Value { Term { [ {V}*30/Pi ] ; In DomainKin ; } } } // Velocity [rpm]
//------------------------------------------------------------------------------------//
//MB2
{ Name P2 ; Value { Term { [ {P2} ] ; In DomainKin2 ; } } } // Position [rad]
{ Name Pdeg2 ; Value { Term { [ {P2}*180/Pi ] ; In DomainKin2 ; } } } // Position [deg]
{ Name V2 ; Value { Term { [ {V2} ] ; In DomainKin2 ; } } } // Velocity [rad/s]
{ Name Vrpm2 ; Value { Term { [ {V2}*30/Pi ] ; In DomainKin2 ; } } } // Velocity [rpm]
}
}
}
```
After that, is just the PostOperation stuff to creat .dat files to plot the results...
I'm sorry for the long issue, but I REALLY can't figure out how to compute those.
Thanks in advance!!!Hello guys! I've been work on a magnetic gear project (https://github.com/CassioKruger/PDD-pure) and I'm having trouble to compute the velocity of the moving bands. I took the "machine_magstadyn_a.pro" file to base my project and added some modifications to use it with 2 rotors instead of 1.
There is a gear ratio between both rotors, so I made this in my .pro file:
```
// pdd gear ratio:
// pH*wH + pL*wL = nP*wP
// wH is the speed of the inner rotor
// wL is the speed of the outer rotor
// wP is the speed of the modulators
// When one of the three parts of the gear is stationary, there will be a constant relation or gear ratio
// between the speeds of other two parts.
// considering that the outer rotor is stationary, the gear ratio becomes:
// -> pH*wH = nP*wP
// -> Gr = pH/nP = wP/wH
// -> gear ratio = nbr of poles at rotor 1 / nbr of modulators
// in this case, the nbr of modulators is equal to the nbr of poles at the outer rotor, so:
gear_ratio = NbrPolesInModel/NbrSectStatorMag;
delta_theta[] = delta_theta_deg * deg2rad ; //rotor 1 step
delta_theta2[] = delta_theta[] * gear_ratio ; //rotor 2 step
```
those "delta_theta[]" should be the step of each rotor to use with "ChangeOfCoordinates" inside the timeloop of the "machine_magstadyn_a_2rotors.pro" file, like this:
```
ChangeOfCoordinates[ NodesOf[Rotor_Moving], RotatePZ[delta_theta[]]];
ChangeOfCoordinates[ NodesOf[Rotor2_Moving], RotatePZ[delta_theta2[]]];
```
and that's working nice, the gear ration is clearly seen at the simulation. But the thing that I'm couldn't figure out is how to compute the correct speed and position of each rotor, so I can show it with graphs and prove that my model is working. I've made this:
-First, I've declared a "DomainKin2", because now there is 2 parts moving
-Then, inside the "Constraint{}" I made this:
```
Constraint{
...
...
...
//Kinetics
{ Name CurrentPosition ; // [m]
Case {
{ Region DomainKin ; Type Init ; Value ($PreviousPosition = 0) ; }
}
}
{ Name CurrentVelocity ; // [rad/s]
Case {
{ Region DomainKin ; Type Init ; Value ($PreviousVelocity = 0) ; }
}
}
//Kinetics - MOVING BAND 2
{ Name CurrentPosition2 ; // [m]
Case {
{ Region DomainKin2 ; Type Init ; Value ($PreviousPosition2 = 0) ; }
}
}
{ Name CurrentVelocity2 ; // [rad/s]
Case {
{ Region DomainKin2 ; Type Init ; Value ($PreviousVelocity2 = 0) ; }
}
}
}
```
-After that, I did the same to "FunctionSpace{}":
```
FunctionSpace{
...
...
...
// For mechanical equation
{ Name Position ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef pr ; Function BF_Region ;
Support DomainKin ; Entity DomainKin ; }
}
GlobalQuantity {
{ Name P ; Type AliasOf ; NameOfCoef pr ; }
}
Constraint {
{ NameOfCoef P ; EntityType Region ; NameOfConstraint CurrentPosition ; }
}
}
{ Name Velocity ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef vr ; Function BF_Region ;
Support DomainKin ; Entity DomainKin ; } }
GlobalQuantity {
{ Name V ; Type AliasOf ; NameOfCoef vr ; }
}
Constraint {
{ NameOfCoef V ; EntityType Region ; NameOfConstraint CurrentVelocity ; }
}
}
// For mechanical equation - MOVING BAND 2
{ Name Position2 ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef pr2 ; Function BF_Region ;
Support DomainKin2 ; Entity DomainKin2 ; }
}
GlobalQuantity {
{ Name P2 ; Type AliasOf ; NameOfCoef pr2 ; }
}
Constraint {
{ NameOfCoef P2 ; EntityType Region ; NameOfConstraint CurrentPosition2 ; }
}
}
{ Name Velocity2 ; Type Scalar ;
BasisFunction {
{ Name sr ; NameOfCoef vr2 ; Function BF_Region ;
Support DomainKin2 ; Entity DomainKin2 ; } }
GlobalQuantity {
{ Name V2 ; Type AliasOf ; NameOfCoef vr2 ; }
}
Constraint {
{ NameOfCoef V2 ; EntityType Region ; NameOfConstraint CurrentVelocity2 ; }
}
}
}
```
-In the "Formulation{}", I did this:
```
Formulation{
...
...
...
// Mechanics
{ Name Mechanical ; Type FemEquation ;
Quantity {
{ Name V ; Type Global ; NameOfSpace Velocity [V] ; } // velocity
{ Name P ; Type Global ; NameOfSpace Position [P] ; } // position
{ Name V2 ; Type Global ; NameOfSpace Velocity2 [V2] ; } // velocity MB2
{ Name P2 ; Type Global ; NameOfSpace Position2 [P2] ; } // position MB2
}
Equation {
GlobalTerm { DtDof [ /*Inertia **/ Dof{V} , {V} ] ; In DomainKin ; }
//GlobalTerm { [ Friction[] * Dof{V} , {V} ] ; In DomainKin ; }
GlobalTerm { [ Torque_mec[], {V} ] ; In DomainKin ; }
GlobalTerm { [ Torque_mag[] , {V} ] ; In DomainKin ; }
GlobalTerm { DtDof [ Dof{P} , {P} ] ; In DomainKin ; }
GlobalTerm { [-Dof{V} , {P} ] ; In DomainKin ; }
//---------------------------------------------------------------//
GlobalTerm { DtDof [ /*Inertia*0.7 **/ Dof{V2} , {V2} ] ; In DomainKin2 ; }
//GlobalTerm { [ Friction[] * Dof{V2} , {V2} ] ; In DomainKin2 ; }
GlobalTerm { [ Torque_mec[] , {V2} ] ; In DomainKin2 ; }
GlobalTerm { [ Torque_mag[] , {V2} ] ; In DomainKin2 ; }
GlobalTerm { DtDof [ Dof{P2} , {P2} ] ; In DomainKin2 ; }
GlobalTerm { [-Dof{V2} , {P2} ] ; In DomainKin2 ; }
}
}
}
```
-And, finally, in the "PostProcessing{}" I did this:
```
PostProcessing{
...
...
...
{ Name Mechanical ; NameOfFormulation Mechanical ;
PostQuantity {
{ Name P ; Value { Term { [ {P} ] ; In DomainKin ; } } } // Position [rad]
{ Name Pdeg ; Value { Term { [ {P}*180/Pi ] ; In DomainKin ; } } } // Position [deg]
{ Name V ; Value { Term { [ {V} ] ; In DomainKin ; } } } // Velocity [rad/s]
{ Name Vrpm ; Value { Term { [ {V}*30/Pi ] ; In DomainKin ; } } } // Velocity [rpm]
//------------------------------------------------------------------------------------//
//MB2
{ Name P2 ; Value { Term { [ {P2} ] ; In DomainKin2 ; } } } // Position [rad]
{ Name Pdeg2 ; Value { Term { [ {P2}*180/Pi ] ; In DomainKin2 ; } } } // Position [deg]
{ Name V2 ; Value { Term { [ {V2} ] ; In DomainKin2 ; } } } // Velocity [rad/s]
{ Name Vrpm2 ; Value { Term { [ {V2}*30/Pi ] ; In DomainKin2 ; } } } // Velocity [rpm]
}
}
}
```
After that, is just the PostOperation stuff to creat .dat files to plot the results...
I'm sorry for the long issue, but I REALLY can't figure out how to compute those.
Thanks in advance!!!https://gitlab.onelab.info/getdp/getdp/issues/49Transient Heat Transfer Problem2018-11-04T21:10:35ZAsimTransient Heat Transfer ProblemHello Everyone,
I am trying to model a transient heat transfer problem (Wellbore within multiple Reservoir Domains of variable properties). I have generated the geometry file and defined the problem. But having some error while running those files (named Cylindrical); like undefined sub-region Reservoir_ID2.
In another case (named Cylindrical 2), I just defined only one reservoir region to overcome this error and it worked. But when I am trying to execute these files, I am getting an error; Not piece-wise Expression: k. for the line where I am defining the thermal conductivity k of my reservoir domain.
I am also trying to generate a GNU plot but it is not generating any file for that and it is only generating the Temp. Map file T.pos.
Can you kindly look into those files and help me resolve these errors. I am a student and a new user of getdp, and doesn't have any coding back-ground, I will appreciate if you can modify the attached files with the resolution of the problem and send them back to me so I can run them.
Best Regards,
Asim Hussain
[Cylindrical.geo](/uploads/0b8c4261e39406160830f5e52e469e6b/Cylindrical.geo)
[Cylindrical.inc](/uploads/da977259ec46d6c5e1d84facca503fec/Cylindrical.inc)
[Cylindrical.pro](/uploads/675dca6aca08591a76ff2fbf498a666a/Cylindrical.pro)
[Cylindrical2.geo](/uploads/73c49a5c8448ae6fd8ba5300a5bc9d19/Cylindrical2.geo)
[Cylindrical2.inc](/uploads/ff58c702c2df2a9c297f08769f8a3366/Cylindrical2.inc)
[Cylindrical2.pro](/uploads/8aea9125042b5ccb3120f78c464cfd04/Cylindrical2.pro)Hello Everyone,
I am trying to model a transient heat transfer problem (Wellbore within multiple Reservoir Domains of variable properties). I have generated the geometry file and defined the problem. But having some error while running those files (named Cylindrical); like undefined sub-region Reservoir_ID2.
In another case (named Cylindrical 2), I just defined only one reservoir region to overcome this error and it worked. But when I am trying to execute these files, I am getting an error; Not piece-wise Expression: k. for the line where I am defining the thermal conductivity k of my reservoir domain.
I am also trying to generate a GNU plot but it is not generating any file for that and it is only generating the Temp. Map file T.pos.
Can you kindly look into those files and help me resolve these errors. I am a student and a new user of getdp, and doesn't have any coding back-ground, I will appreciate if you can modify the attached files with the resolution of the problem and send them back to me so I can run them.
Best Regards,
Asim Hussain
[Cylindrical.geo](/uploads/0b8c4261e39406160830f5e52e469e6b/Cylindrical.geo)
[Cylindrical.inc](/uploads/da977259ec46d6c5e1d84facca503fec/Cylindrical.inc)
[Cylindrical.pro](/uploads/675dca6aca08591a76ff2fbf498a666a/Cylindrical.pro)
[Cylindrical2.geo](/uploads/73c49a5c8448ae6fd8ba5300a5bc9d19/Cylindrical2.geo)
[Cylindrical2.inc](/uploads/ff58c702c2df2a9c297f08769f8a3366/Cylindrical2.inc)
[Cylindrical2.pro](/uploads/8aea9125042b5ccb3120f78c464cfd04/Cylindrical2.pro)https://gitlab.onelab.info/getdp/getdp/issues/48Wrench2D Tutorial: Second order elements & inhomogeneous Dirichlet data2018-10-24T08:45:08ZMarc Alexander SchweitzerWrench2D Tutorial: Second order elements & inhomogeneous Dirichlet dataI made some minor changes to the wrench2D tutorial to have inhomogeneous Dirichlet data (prescribed non-zero displacement). It seems to work for linear elements but when activating the second order option/flag, the results near the Dirichlet boundary with inhomogeneous data are wrong. I attached the updated wrench2D.pro as well as screen shots of the results with linear / second order. Are further changes required to deal with inhomogeneous data? Or is there an issue with the elimination of essential boundary data when using higher order elements.!
[linear_elements](/uploads/b609dbf1d29882eb3bf4b0b21da6fdb8/linear_elements.png)
![second_order_elements](/uploads/b3800fe006f0e949abfba82efd1079e3/second_order_elements.png)
[wrench2D.pro](/uploads/52a50262d45c4be27ec6287d3f755627/wrench2D.pro)I made some minor changes to the wrench2D tutorial to have inhomogeneous Dirichlet data (prescribed non-zero displacement). It seems to work for linear elements but when activating the second order option/flag, the results near the Dirichlet boundary with inhomogeneous data are wrong. I attached the updated wrench2D.pro as well as screen shots of the results with linear / second order. Are further changes required to deal with inhomogeneous data? Or is there an issue with the elimination of essential boundary data when using higher order elements.!
[linear_elements](/uploads/b609dbf1d29882eb3bf4b0b21da6fdb8/linear_elements.png)
![second_order_elements](/uploads/b3800fe006f0e949abfba82efd1079e3/second_order_elements.png)
[wrench2D.pro](/uploads/52a50262d45c4be27ec6287d3f755627/wrench2D.pro)https://gitlab.onelab.info/getdp/getdp/issues/42Inductor model2018-05-04T10:42:08ZMattia ConteInductor modelHi there,
I am trying to use your "Inductor" model for obtaining the inductance in Time Domain simulation.
First, in the formula to calculate the inductance you should substitute the term "II" with the actual value of the current "$I"
to obtain the good result.
However, the simulation gives wrong results (in 2D and 3D) simulating the system with a conducting core.
The inductance in function of time should be constant, and should be lower with a conducting core because eddy currents creates
an opposite flux that make the total flux (then the inductance L=Psi/I) decrease.
Do you know what can be the cause of it ?
Thanks in advance,
MattiaHi there,
I am trying to use your "Inductor" model for obtaining the inductance in Time Domain simulation.
First, in the formula to calculate the inductance you should substitute the term "II" with the actual value of the current "$I"
to obtain the good result.
However, the simulation gives wrong results (in 2D and 3D) simulating the system with a conducting core.
The inductance in function of time should be constant, and should be lower with a conducting core because eddy currents creates
an opposite flux that make the total flux (then the inductance L=Psi/I) decrease.
Do you know what can be the cause of it ?
Thanks in advance,
Mattiahttps://gitlab.onelab.info/getdp/getdp/issues/38Segmentation fault in GetDP 2.11.0 64-bit Windows2017-08-05T11:55:13ZChristophe GeuzaineSegmentation fault in GetDP 2.11.0 64-bit WindowsWorking with the enclosed files, the command _getdp 2D_eucard_v3.pro -pre ResolutionH -msh 2D_eucard_v3.msh_ behaves differently between the 32 and 64-bit versions of the GetDP 2.11.0 Windows executable. With the 32-bit version, the command works fine whereas with the 64-bit executable it generates a segmentation fault.
The error seems to happen while executing
```
GlobalQuantity_P = (struct GlobalQuantity*)
List_Pointer(QuantityStorage_P->FunctionSpace->GlobalQuantity,
*(int*) List_Pointer(DefineQuantity_P->IndexInFunctionSpace,0)
```
in file Kernel/Treatment_Formulation.cpp around line 669. In the 64-bit Windows version, the pointer to QuantityStorage_P->FunctionSpace is invalid (0xFFFFFFFF00000001) while it was ok in the 32-bit version. I figured the problem could be circumvented by increasing NBR_MAX_BASISFUNCTIONS in Interface/ProData.h.
[2D_eucard_v3.pro](/uploads/78e938cdc91e79d86c6b603a087ce45b/2D_eucard_v3.pro)
[2D_eucard_parameters_v3.geo](/uploads/16420bfc22441f22d67e815eafed88d2/2D_eucard_parameters_v3.geo)
[2D_eucard_v3.geo](/uploads/9f71660fedf134898d5127cc0c20e3ce/2D_eucard_v3.geo)
[2D_eucard_macros_v3.geo](/uploads/b505c587edded44d6cb5ea5e40f2192e/2D_eucard_macros_v3.geo)
[2D_eucard_display_v3.geo](/uploads/3713f5ceb202387d416462033a4a51ca/2D_eucard_display_v3.geo)
Working with the enclosed files, the command _getdp 2D_eucard_v3.pro -pre ResolutionH -msh 2D_eucard_v3.msh_ behaves differently between the 32 and 64-bit versions of the GetDP 2.11.0 Windows executable. With the 32-bit version, the command works fine whereas with the 64-bit executable it generates a segmentation fault.
The error seems to happen while executing
```
GlobalQuantity_P = (struct GlobalQuantity*)
List_Pointer(QuantityStorage_P->FunctionSpace->GlobalQuantity,
*(int*) List_Pointer(DefineQuantity_P->IndexInFunctionSpace,0)
```
in file Kernel/Treatment_Formulation.cpp around line 669. In the 64-bit Windows version, the pointer to QuantityStorage_P->FunctionSpace is invalid (0xFFFFFFFF00000001) while it was ok in the 32-bit version. I figured the problem could be circumvented by increasing NBR_MAX_BASISFUNCTIONS in Interface/ProData.h.
[2D_eucard_v3.pro](/uploads/78e938cdc91e79d86c6b603a087ce45b/2D_eucard_v3.pro)
[2D_eucard_parameters_v3.geo](/uploads/16420bfc22441f22d67e815eafed88d2/2D_eucard_parameters_v3.geo)
[2D_eucard_v3.geo](/uploads/9f71660fedf134898d5127cc0c20e3ce/2D_eucard_v3.geo)
[2D_eucard_macros_v3.geo](/uploads/b505c587edded44d6cb5ea5e40f2192e/2D_eucard_macros_v3.geo)
[2D_eucard_display_v3.geo](/uploads/3713f5ceb202387d416462033a4a51ca/2D_eucard_display_v3.geo)
https://gitlab.onelab.info/getdp/getdp/issues/37magnetostatics simulation gives wrong results2017-08-05T11:56:16ZChristophe Geuzainemagnetostatics simulation gives wrong resultsI made a model with one magnet and an iron frame (mu_r = 9000) wrapped around it. I expect most of the B lines should go inside the iron frame, but the simulation result seems as if the iron frame does not exist. Please see attached from the geo model. Open the geo model in gmsh first, then merge (gmsh menu: File >> Merge) in magnetostatics.pro in the template folder came with gmsh. You should be able to set the model interactively. Set a constant magnetization for the magnet of 900000 in z direction. Select the right materials for air, set frame with constant mu_r = 9000, and the boundary condition on "Inf" can be either way (may leave at its default value). Then click Run.
Below is the 1magnet.geo file I used:
```
//=================start==================
// define geometry-specific parameters
mm = 1.e-3;
DefineConstant[
cub = {10*mm, Name "Parameters/2Magnet bottom size [m]"}
hite = {20*mm, Name "Parameters/2Magnet hieght [m]"}
lc1 = {TotalMemory <= 2048 ? 5*mm : 2*mm, Name "Parameters/3Mesh size on magnets [m]"}
lc2 = {TotalMemory <= 2048 ? 20*mm : 10*mm, Name "Parameters/4Mesh size at infinity [m]"}
inf = {100*mm, Name "Parameters/1Air box distance [m]"}
];
// change global Gmsh options
Mesh.Optimize = 1; // optimize quality of tetrahedra
Mesh.VolumeEdges = 0; // hide volume edges
Geometry.ExactExtrusion = 0; // to allow rotation of extruded shapes
Solver.AutoMesh = 2; // always remesh if necessary (don't reuse mesh on disk)
p1 = newp; Point(p1) = {-cub, -cub, -hite, lc1};
p2 = newp; Point(p2) = { cub, -cub, -hite, lc1};
p3 = newp; Point(p3) = { cub, cub, -hite, lc1};
p4 = newp; Point(p4) = {-cub, cub, -hite, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};
s1 = news; Plane Surface(s1) = {ll1};
mag[] = Extrude {0, 0, 2*hite} { Surface{s1}; };
Physical Volume("Magnet") = {mag[1]};
//create steel frame around the magnet
p1 = newp; Point(p1) = {-2*cub, -cub, -hite, lc1};
p2 = newp; Point(p2) = { 2*cub, -cub, -hite, lc1};
p3 = newp; Point(p3) = { 2*cub, -cub, hite, lc1};
p4 = newp; Point(p4) = {-2*cub, -cub, hite, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};
hite2 = hite + cub;
p1 = newp; Point(p1) = {-4*cub, -cub, -hite2, lc1};
p2 = newp; Point(p2) = { 4*cub, -cub, -hite2, lc1};
p3 = newp; Point(p3) = { 4*cub, -cub, hite2, lc1};
p4 = newp; Point(p4) = {-4*cub, -cub, hite2, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll2 = newll; Line Loop(ll2) = {l1,l2,l3,l4};
s1 = news; Plane Surface(s1) = {ll2, ll1};
frame[] = Extrude {0, 2*cub, 0} { Surface{s1}; };
Physical Volume("Frame") = {frame[1]};
// create air box around magnets
BoundingBox; // recompute model bounding box
cx = (General.MinX + General.MaxX) / 2;
cy = (General.MinY + General.MaxY) / 2;
cz = (General.MinZ + General.MaxZ) / 2;
lx = 2*inf + General.MaxX - General.MinX;
ly = 2*inf + General.MaxY - General.MinZ;
lz = 2*inf + General.MaxZ - General.MinZ;
p1 = newp; Point (p1) = {cx-lx/2, cy-ly/2, cz-lz/2, lc2};
p2 = newp; Point (p2) = {cx+lx/2, cy-ly/2, cz-lz/2, lc2};
l1 = newl; Line(l1) = {p1, p2};
e1[] = Extrude {0, ly, 0} { Line{l1}; };
air[] = Extrude {0, 0, lz} { Surface{e1[1]}; };
slair = newsl; Surface Loop(slair) = Boundary{Volume{air[1]}; };
slmag = newsl; Surface Loop(slmag) = Boundary{Volume{mag[1]}; };
vair = newv; Volume(vair) = {slair, slmag};
Physical Volume("Air") = vair; // air
Physical Surface("Inf") = {e1[1], air[0], air[2], air[3], air[4], air[5]};
Delete { Volume{air[1]}; }
//=================end==================
```
[1magnet.geo](/uploads/7f08801b5b3492eb6d959a7606e7d989/1magnet.geo)I made a model with one magnet and an iron frame (mu_r = 9000) wrapped around it. I expect most of the B lines should go inside the iron frame, but the simulation result seems as if the iron frame does not exist. Please see attached from the geo model. Open the geo model in gmsh first, then merge (gmsh menu: File >> Merge) in magnetostatics.pro in the template folder came with gmsh. You should be able to set the model interactively. Set a constant magnetization for the magnet of 900000 in z direction. Select the right materials for air, set frame with constant mu_r = 9000, and the boundary condition on "Inf" can be either way (may leave at its default value). Then click Run.
Below is the 1magnet.geo file I used:
```
//=================start==================
// define geometry-specific parameters
mm = 1.e-3;
DefineConstant[
cub = {10*mm, Name "Parameters/2Magnet bottom size [m]"}
hite = {20*mm, Name "Parameters/2Magnet hieght [m]"}
lc1 = {TotalMemory <= 2048 ? 5*mm : 2*mm, Name "Parameters/3Mesh size on magnets [m]"}
lc2 = {TotalMemory <= 2048 ? 20*mm : 10*mm, Name "Parameters/4Mesh size at infinity [m]"}
inf = {100*mm, Name "Parameters/1Air box distance [m]"}
];
// change global Gmsh options
Mesh.Optimize = 1; // optimize quality of tetrahedra
Mesh.VolumeEdges = 0; // hide volume edges
Geometry.ExactExtrusion = 0; // to allow rotation of extruded shapes
Solver.AutoMesh = 2; // always remesh if necessary (don't reuse mesh on disk)
p1 = newp; Point(p1) = {-cub, -cub, -hite, lc1};
p2 = newp; Point(p2) = { cub, -cub, -hite, lc1};
p3 = newp; Point(p3) = { cub, cub, -hite, lc1};
p4 = newp; Point(p4) = {-cub, cub, -hite, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};
s1 = news; Plane Surface(s1) = {ll1};
mag[] = Extrude {0, 0, 2*hite} { Surface{s1}; };
Physical Volume("Magnet") = {mag[1]};
//create steel frame around the magnet
p1 = newp; Point(p1) = {-2*cub, -cub, -hite, lc1};
p2 = newp; Point(p2) = { 2*cub, -cub, -hite, lc1};
p3 = newp; Point(p3) = { 2*cub, -cub, hite, lc1};
p4 = newp; Point(p4) = {-2*cub, -cub, hite, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};
hite2 = hite + cub;
p1 = newp; Point(p1) = {-4*cub, -cub, -hite2, lc1};
p2 = newp; Point(p2) = { 4*cub, -cub, -hite2, lc1};
p3 = newp; Point(p3) = { 4*cub, -cub, hite2, lc1};
p4 = newp; Point(p4) = {-4*cub, -cub, hite2, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
ll2 = newll; Line Loop(ll2) = {l1,l2,l3,l4};
s1 = news; Plane Surface(s1) = {ll2, ll1};
frame[] = Extrude {0, 2*cub, 0} { Surface{s1}; };
Physical Volume("Frame") = {frame[1]};
// create air box around magnets
BoundingBox; // recompute model bounding box
cx = (General.MinX + General.MaxX) / 2;
cy = (General.MinY + General.MaxY) / 2;
cz = (General.MinZ + General.MaxZ) / 2;
lx = 2*inf + General.MaxX - General.MinX;
ly = 2*inf + General.MaxY - General.MinZ;
lz = 2*inf + General.MaxZ - General.MinZ;
p1 = newp; Point (p1) = {cx-lx/2, cy-ly/2, cz-lz/2, lc2};
p2 = newp; Point (p2) = {cx+lx/2, cy-ly/2, cz-lz/2, lc2};
l1 = newl; Line(l1) = {p1, p2};
e1[] = Extrude {0, ly, 0} { Line{l1}; };
air[] = Extrude {0, 0, lz} { Surface{e1[1]}; };
slair = newsl; Surface Loop(slair) = Boundary{Volume{air[1]}; };
slmag = newsl; Surface Loop(slmag) = Boundary{Volume{mag[1]}; };
vair = newv; Volume(vair) = {slair, slmag};
Physical Volume("Air") = vair; // air
Physical Surface("Inf") = {e1[1], air[0], air[2], air[3], air[4], air[5]};
Delete { Volume{air[1]}; }
//=================end==================
```
[1magnet.geo](/uploads/7f08801b5b3492eb6d959a7606e7d989/1magnet.geo)https://gitlab.onelab.info/getdp/getdp/issues/30OnGrid Interpolation with Timestepping results in erronous output for the fir...2017-03-27T09:20:47ZChristophe GeuzaineOnGrid Interpolation with Timestepping results in erronous output for the first tilmestepIf I use OnGrid interpolation with TimeStepping in PostOperation, the first timestep that is outputted has some small errors. I've modified the Simple_RLC example to demonstrate this.
I output timesteps 90:101 and 91:101 of the vector current. When subtracting two equal time steps the difference can be seen for any subtraction that includes the first time step. Subtractions of later time steps results in zero as expected. See the attached screenshot.If I use OnGrid interpolation with TimeStepping in PostOperation, the first timestep that is outputted has some small errors. I've modified the Simple_RLC example to demonstrate this.
I output timesteps 90:101 and 91:101 of the vector current. When subtracting two equal time steps the difference can be seen for any subtraction that includes the first time step. Subtractions of later time steps results in zero as expected. See the attached screenshot.https://gitlab.onelab.info/getdp/getdp/issues/27time-dependent functions with time stepping schemes other than implicit Euler2017-03-27T09:20:45ZChristophe Geuzainetime-dependent functions with time stepping schemes other than implicit EulerGetDP automatically handles time-dependent constraints when they are provided using the TimeFunction mechanism in an Assign-type Constraint.
However, GetDP cannot automatically transform general time-dependent source terms in weak formulations. Such source terms will be correctly treated only for implicit Euler, as the expression in the Galerkin term is evaluated at the current time step.
For other schemes, the source term should be written explicitly, by splitting it in two (theta f_n+1 + (1-theta) f_n), making use of the AtAnteriorTimeStep[] for the second part, and specifying NeverDt in the Galerkin term.GetDP automatically handles time-dependent constraints when they are provided using the TimeFunction mechanism in an Assign-type Constraint.
However, GetDP cannot automatically transform general time-dependent source terms in weak formulations. Such source terms will be correctly treated only for implicit Euler, as the expression in the Galerkin term is evaluated at the current time step.
For other schemes, the source term should be written explicitly, by splitting it in two (theta f_n+1 + (1-theta) f_n), making use of the AtAnteriorTimeStep[] for the second part, and specifying NeverDt in the Galerkin term.https://gitlab.onelab.info/getdp/getdp/issues/18Add function that gives size of geometric bounding box2017-03-27T09:20:40ZChristophe GeuzaineAdd function that gives size of geometric bounding boxWe should add a function that returns the values of the geom bbox stored in
GeoData_P->Xmin, etc.
This would be useful for determining regularization constantsWe should add a function that returns the values of the geom bbox stored in
GeoData_P->Xmin, etc.
This would be useful for determining regularization constantshttps://gitlab.onelab.info/getdp/getdp/issues/10generalize localterm2017-08-05T12:01:19ZChristophe Geuzainegeneralize localtermgeneralize localterm (equation part should call Cal_vBFxDof)generalize localterm (equation part should call Cal_vBFxDof)https://gitlab.onelab.info/getdp/getdp/issues/9Should recompute Current.x,y,z in Cal_vBFxDof?2017-03-27T09:20:54ZChristophe GeuzaineShould recompute Current.x,y,z in Cal_vBFxDof?subject says it all...subject says it all...https://gitlab.onelab.info/getdp/getdp/issues/6SetFrequency2017-03-27T09:20:53ZChristophe GeuzaineSetFrequencyverify if the post-pro uses the correct frequency for each system if
the frequency was changed by hand using SetFrequency in the
resolution, in between two GenerateSystemverify if the post-pro uses the correct frequency for each system if
the frequency was changed by hand using SetFrequency in the
resolution, in between two GenerateSystemhttps://gitlab.onelab.info/getdp/getdp/issues/4access Time & TimeImag in post-processing2017-03-27T09:20:53ZChristophe Geuzaineaccess Time & TimeImag in post-processingSubject says it all (mostly useful for eigenvalues Re/Im)Subject says it all (mostly useful for eigenvalues Re/Im)https://gitlab.onelab.info/getdp/getdp/issues/3complex conjugation problem with sparskit2017-03-27T09:20:53ZChristophe Geuzainecomplex conjugation problem with sparskitWith sparsekit and 3D modelling in getdp I found an
error. If one includes terms like
Galerkin { [ Einc[], {E} ]; In port; Integration I1;
Jacobian Jac;}
Where Einc is purely imaginary then sparsekit complex conjugates it.
It is only a problem in 3D. In 2D it does not exist. The problem
also disappers with petsc.With sparsekit and 3D modelling in getdp I found an
error. If one includes terms like
Galerkin { [ Einc[], {E} ]; In port; Integration I1;
Jacobian Jac;}
Where Einc is purely imaginary then sparsekit complex conjugates it.
It is only a problem in 3D. In 2D it does not exist. The problem
also disappers with petsc.