Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

How do Variable constraints Work? #532

Open
klbudzin opened this issue Feb 1, 2018 · 14 comments
Open

How do Variable constraints Work? #532

klbudzin opened this issue Feb 1, 2018 · 14 comments

Comments

@klbudzin
Copy link
Collaborator

klbudzin commented Feb 1, 2018

I am currently trying to pin down my Temperature to 500 at an x position of 0.05 in my One dimensional example. I'm not quite sure how constraining points work, and attempted setting a constraint point but ran into an error.

My variables section in my input file is:

[Variables]
[./Temperature]
names = 'T'
fe_family = 'LAGRANGE'
order = 'FIRST'

  [./ConstrainedPoints]
[./FirstEdgeNode]
  constraint_locations = '0.05'
  constraining_points_x = '0.05'
  constraining_points_coeff = '1'
  constraining_points_var = 'T'
  constraining_rhs = '500'
    [../]
  [../]		   

[../SingleVariable]
names = 'M_dot'
fe_family = 'LAGRANGE'
order = 'FIRST'
[../SpeciesMassFractions]
names = 'Y_'
fe_family = 'LAGRANGE'
order = 'FIRST'
material = 'OzoneGas'
[]

The error I was obtaining was:

ailed to find a Node at point (x,y,z)=( 0, 0, 0)
[Thread debugging using libthread_db enabled]
0x000000311deac82e in waitpid () from /lib64/libc.so.6
To enable execution of this file add
add-auto-load-safe-path /zoidberg1/data/shared/software/apps/gcc/7.2.0/lib64/libstdc++.so.6.0.24-gdb.py
line to your configuration file "/nsm/home/klbudzin/.gdbinit".
To completely disable this security protection add
set auto-load safe-path /
line to your configuration file "/nsm/home/klbudzin/.gdbinit".
For more information about this security protection see the
"Auto-loading safe path" section in the GDB manual. E.g., run from the shell:
info "(gdb)Auto-loading safe path"
#0 0x000000311deac82e in waitpid () from /lib64/libc.so.6
#1 0x000000311de3e479 in do_system () from /lib64/libc.so.6
#2 0x00007f662a9506a4 in libMesh::print_trace(std::basic_ostream<char, std::char_traits >&) () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#3 0x00007f662a94cbf7 in libMesh::MacroFunctions::report_error(char const*, int, char const*, char const*) () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#4 0x00007f662bae5e91 in GRINS::ConstrainedPoints::constrain() () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#5 0x00007f662b0a0448 in libMesh::System::reinit_constraints() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#6 0x00007f662b0a7eec in libMesh::System::init_data() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#7 0x00007f662b088e54 in libMesh::ImplicitSystem::init_data() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#8 0x00007f662b0596e5 in libMesh::DifferentiableSystem::init_data() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#9 0x00007f662baf07dd in GRINS::MultiphysicsSystem::init_data() () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#10 0x00007f662b0a81a6 in libMesh::System::init() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#11 0x00007f662b06216d in libMesh::EquationSystems::init() () from /zoidberg1/data/shared/software/libs/libmesh/master/gcc/7.2.0/mpich/3.2.1/petsc/3.8.3/openblas/0.2.20/boost/1.66.0/hdf5/1.10.1/vtk/8.1.0/lib/libmesh_opt.so.0
#12 0x00007f662be8246f in GRINS::Solver::initialize(GetPot const&, std::shared_ptrlibMesh::EquationSystems, GRINS::MultiphysicsSystem*) () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#13 0x00007f662be904f9 in GRINS::Simulation::init_multiphysics_system(GetPot const&) () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#14 0x00007f662be93de6 in GRINS::Simulation::Simulation(GetPot const&, GetPot&, GRINS::SimulationBuilder&, libMesh::Parallel::Communicator const&) () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#15 0x00007f662beca262 in GRINS::Runner::init() () from /zoidberg1/data/shared/klbudzin/WorkingGrins/OptBuild/lib/libgrins.so.0
#16 0x000000000040297d in main ()
[0] ../../grins/src/constraints/src/constrained_points.C, line 205, compiled Jan 31 2018 at 12:49:56

Thanks, @roystgnr

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

Hmm... clearly I needed to document the ConstrainedPoints class better. Or at all.

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

But the first problem here appears to be that you have a constraint_locations input options, whereas ConstrainedPoints will be looking for constraint_location.

So you're trying to pin a node at 0.05,0,0, but the ConstrainedPoints code is looking for and failing to find a node at the default 0,0,0 location.

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

You're also setting up an inconsistent constraint: with the constraining_points location equal to the constrained point, and a coefficient of 1, and the same variable, and an rhs of 500, you're trying to impose T(.05) = 1*T(.05)+500. This will not happen. If you're lucky then you'll run a new git libMesh in debug or devel mode, and the new constraint sanity checking @dknez added to libMesh in libMesh/libmesh#1562 will scream and die before silently handing back incorrect results.

@pbauman
Copy link
Member

pbauman commented Feb 1, 2018

But the first problem here appears to be that you have a constraint_locations input options, whereas ConstrainedPoints will be looking for constraint_location.

Can we set this to be a required option, i.e. error with message constraint_location not found or is it an optional parameter? If required, @klbudzin or I can drop in a quick PR.

@pbauman
Copy link
Member

pbauman commented Feb 1, 2018

You're also setting up an inconsistent constraint: with the constraining_points location equal to the constrained point, and a coefficient of 1, and the same variable, and an rhs of 500, you're trying to impose T(.05) = 1*T(.05)+500.

Clearly we're not understanding how the constraint code works. Can you quick us a quick rundown here so we can push a doc PR or could you push a quick doc PR?

@pbauman
Copy link
Member

pbauman commented Feb 1, 2018

Clearly we're not understanding
We = @klbudzin and I

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

@pbauman, do we actually have a GRINS capability to enforce a single variable constraint with a hard constraint equation?

Based on your replies, the answer is "no"? Then the closest we can come right now would be the VariablePinning class, which does pinning via a penalty method. We usually use it to pin pressure values, but aside from the atavistic name in the PressurePinning helper class, there's nothing left that's pressure-specific.

I don't see any examples of it in our input files, but I used it all the time for flame experiments, e.g.:

[Physics]
   [./VariablePinning]
     penalty = 1e4
     pin_variable = 'T'
     pin_value = 1000
     pin_location = '0 0 0'
   [../] 
[]     

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

constraint_location shouldn't be an optional parameter; the utility of allowing someone to default to 0,0,0 on purpose with one less line of config file is almost certainly outweighted by the disutility of not giving a sane error message when someone defaults to 0,0,0 by accident.

@pbauman
Copy link
Member

pbauman commented Feb 1, 2018

do we actually have a GRINS capability to enforce a single variable constraint with a hard constraint equation?

I guess I thought that's what the ConstrainedPoints was. I'm inferring from your statement that in fact it's not that. :)

@pbauman
Copy link
Member

pbauman commented Feb 1, 2018

constraint_location shouldn't be an optional parameter;

OK, we'll make a quick fix.

the utility of allowing someone to default to 0,0,0 on purpose with one less line of config file is almost certainly outweighted by the disutility of not giving a sane error message when someone defaults to 0,0,0 by accident.

Yes. I've learned this lesson too many times now. :)

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

The closest thing I have to documentation is test/input_files/poisson_user_constraint.in

Basically:

Each ConstrainedPoints subsection is a new point whose value for that variable will be constrained.

constraint_location is x,y,z coordinates for the constrained point.
constraining_points_x (y,z) are vectors of coordinates of points constraining it.
constraining_points_coeff is a vector of coefficients in the constraint equation
constraining_points_var is a vector of which variables are used on each constraining point
constraint_rhs is a constant offset.

Huh. I guess you actually could use ConstrainedPoints for what you want, just by leaving all those vectors empty! I've never tried it so there might be some sanity checking which prevents it from working but if so that'd be easy to fix.

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

Oh, and forbid_overwrite=false lets you override some pre-existing (e.g. Dirichlet or periodic boundary) constraint equation. If you have to use this then it's strong evidence that something is wrong with your mathematical model, but the software can be told that you're the boss.

@roystgnr
Copy link
Member

roystgnr commented Feb 1, 2018

I don't have time to write a proper docs PR at the moment, but I'll review/edit anything you guys want to put together if you're in a hurry.

@pbauman
Copy link
Member

pbauman commented Mar 1, 2018

Just FYI, this is all working as expected, I'm just keeping open to submit a doc PR.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants