getdp issues https://gitlab.onelab.info/getdp/getdp/-/issues 2019-08-01T05:41:14Z https://gitlab.onelab.info/getdp/getdp/-/issues/57 Link Constraint for Non-Identical Regions 2019-08-01T05:41:14Z Bruno de Sousa Alves Link Constraint for Non-Identical Regions Dear 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 Alves Dear 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 Alves https://gitlab.onelab.info/getdp/getdp/-/issues/55 Computing rotor speed and position 2019-08-01T05:44:30Z Cássio Kruger krugercassio@gmail.com Computing rotor speed and position 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!!! 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/49 Transient Heat Transfer Problem 2018-11-04T21:10:35Z Asim Transient Heat Transfer Problem 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) 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/48 Wrench2D Tutorial: Second order elements &amp; inhomogeneous Dirichlet data 2018-10-24T08:45:08Z Marc Alexander Schweitzer Wrench2D Tutorial: Second order elements & inhomogeneous Dirichlet data 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) 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/42 Inductor model 2018-05-04T10:42:08Z Mattia Conte Inductor model Hi 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, Mattia Hi 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, Mattia https://gitlab.onelab.info/getdp/getdp/-/issues/38 Segmentation fault in GetDP 2.11.0 64-bit Windows 2017-08-05T11:55:13Z Christophe Geuzaine Segmentation fault in GetDP 2.11.0 64-bit Windows 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) 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/37 magnetostatics simulation gives wrong results 2017-08-05T11:56:16Z Christophe Geuzaine magnetostatics simulation gives wrong results 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}; //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}; // 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}; }; slair = newsl; Surface Loop(slair) = Boundary{Volume{air}; }; slmag = newsl; Surface Loop(slmag) = Boundary{Volume{mag}; }; vair = newv; Volume(vair) = {slair, slmag}; Physical Volume("Air") = vair; // air Physical Surface("Inf") = {e1, air, air, air, air, air}; Delete { Volume{air}; } //=================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}; //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}; // 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}; }; slair = newsl; Surface Loop(slair) = Boundary{Volume{air}; }; slmag = newsl; Surface Loop(slmag) = Boundary{Volume{mag}; }; vair = newv; Volume(vair) = {slair, slmag}; Physical Volume("Air") = vair; // air Physical Surface("Inf") = {e1, air, air, air, air, air}; Delete { Volume{air}; } //=================end================== ``` [1magnet.geo](/uploads/7f08801b5b3492eb6d959a7606e7d989/1magnet.geo) https://gitlab.onelab.info/getdp/getdp/-/issues/30 OnGrid Interpolation with Timestepping results in erronous output for the fir... 2017-03-27T09:20:47Z Christophe Geuzaine OnGrid Interpolation with Timestepping results in erronous output for the first tilmestep 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. 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/27 time-dependent functions with time stepping schemes other than implicit Euler 2017-03-27T09:20:45Z Christophe Geuzaine time-dependent functions with time stepping schemes other than implicit Euler 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. 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/18 Add function that gives size of geometric bounding box 2017-03-27T09:20:40Z Christophe Geuzaine Add function that gives size of geometric bounding box We 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 constants We 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 constants https://gitlab.onelab.info/getdp/getdp/-/issues/10 generalize localterm 2017-08-05T12:01:19Z Christophe Geuzaine generalize localterm generalize localterm (equation part should call Cal_vBFxDof) generalize localterm (equation part should call Cal_vBFxDof) https://gitlab.onelab.info/getdp/getdp/-/issues/9 Should recompute Current.x,y,z in Cal_vBFxDof? 2017-03-27T09:20:54Z Christophe Geuzaine Should recompute Current.x,y,z in Cal_vBFxDof? subject says it all... subject says it all... https://gitlab.onelab.info/getdp/getdp/-/issues/6 SetFrequency 2021-02-21T16:36:38Z Christophe Geuzaine SetFrequency verify 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 GenerateSystem verify 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 GenerateSystem https://gitlab.onelab.info/getdp/getdp/-/issues/4 access Time &amp; TimeImag in post-processing 2017-03-27T09:20:53Z Christophe Geuzaine access Time & TimeImag in post-processing Subject 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/3 complex conjugation problem with sparskit 2020-05-29T15:44:21Z Christophe Geuzaine complex conjugation problem with sparskit 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. 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.