#398 not reproducible on Gridap v0.14.0

[bhaveshshrimali]

Hi @fverdugo @santiagobadia

I just installed Julia 1.5.1 and Gridap 0.14.0 on my cluster and tried running the example code from #398 . It is failing with

julia dolfinHyperelasticity.jl
ERROR: LoadError: UndefVarError: Args not defined
Stacktrace:
 [1] Base.Broadcast.Broadcasted{Base.Broadcast.Style{Tuple},Axes,F,Args} where Args<:Tuple where F where Axes(::typeof(==), ::Tuple{Tuple{var"#s313"} where var"#s313"<:Int64,Int64}, ::Nothing) at ./broadcast.jl:179 (repeats 2 times)
 [2] broadcasted at ./broadcast.jl:1265 [inlined]
 [3] broadcasted at ./broadcast.jl:1263 [inlined]
 [4] _admissible_permutations(::Gridap.ReferenceFEs.DFace{2}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/ExtrusionPolytopes.jl:683
 [5] _precompute_vertex_perms_if_possible(::Gridap.ReferenceFEs.DFace{2}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/ExtrusionPolytopes.jl:670
 [6] Gridap.ReferenceFEs.ExtrusionPolytope(::Gridap.ReferenceFEs.DFace{2}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/ExtrusionPolytopes.jl:53
 [7] Polytope at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/ExtrusionPolytopes.jl:146 [inlined]
 [8] #3 at ./none:0 [inlined]
 [9] iterate at ./generator.jl:47 [inlined]
 [10] collect(::Base.Generator{UnitRange{Int64},Gridap.ReferenceFEs.var"#3#4"{2,Gridap.ReferenceFEs.ExtrusionPolytope{3}}}) at ./array.jl:686
 [11] _compute_reffaces_and_face_types at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/Polytopes.jl:592 [inlined]
 [12] get_reffaces(::Type{Polytope{2}}, ::Gridap.ReferenceFEs.ExtrusionPolytope{3}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/ReferenceFEs/Polytopes.jl:552
 [13] #7 at ./none:0 [inlined]
 [14] iterate at ./generator.jl:47 [inlined]
 [15] collect(::Base.Generator{Array{Gridap.ReferenceFEs.ExtrusionPolytope{3},1},Gridap.Geometry.var"#7#9"{2}}) at ./array.jl:686
 [16] compute_reffaces(::Type{Polytope{2}}, ::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/GridTopologies.jl:318
 [17] get_face_type(::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}, ::Int64) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/GridTopologies.jl:308
 [18] _setup_nface_to_mface!(::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}, ::Int64, ::Int64) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/UnstructuredGridTopologies.jl:309
 [19] _setup_faces!(::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}, ::Int64, ::Int64) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/UnstructuredGridTopologies.jl:199
 [20] get_faces at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/UnstructuredGridTopologies.jl:164 [inlined]
 [21] _fix_geolabels(::Int64, ::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}, ::Array{Array{Int32,1},1}, ::Int64, ::Int64) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:137
 [22] _fill_cartesian_entities!(::Gridap.Geometry.FaceLabeling, ::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:130
 [23] _fill_cartesian_face_labeling!(::Gridap.Geometry.FaceLabeling, ::Gridap.Geometry.UnstructuredGridTopology{3,3,Float64,true}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:108
 [24] CartesianDiscreteModel(::Gridap.Geometry.CartesianDescriptor{3,Float64,typeof(identity)}) at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:26
 [25] #CartesianDiscreteModel#84 at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:69 [inlined]
 [26] CartesianDiscreteModel at /home/bshrima2/.julia/packages/Gridap/y0NoW/src/Geometry/CartesianDiscreteModels.jl:68 [inlined]
 [27] run() at /projects/meca/bshrima2/scikitFEMCodes/dolfinHyperelasticity.jl:50
 [28] top-level scope at /projects/meca/bshrima2/scikitFEMCodes/dolfinHyperelasticity.jl:100
 [29] include(::Function, ::Module, ::String) at ./Base.jl:380
 [30] include(::Module, ::String) at ./Base.jl:368
 [31] exec_options(::Base.JLOptions) at ./client.jl:296
 [32] _start() at ./client.jl:506
in expression starting at /projects/meca/bshrima2/scikitFEMCodes/dolfinHyperelasticity.jl:100

Downgrading back to v0.13.4 which I had been using on my laptop seems to run without any problems. Do you think it is an issue with my environment specifically or something else?

Here is the code just as reference (it's basically the same as in #398 but with simplexify)

Code

using Gridap
using LinearAlgebra: inv, det
using LineSearches: BackTracking, StrongWolfe, HagerZhang
using BenchmarkTools

E = 10.
ν = 0.3
const μ = E/2/(1+ν)
const λ = 2*μ*ν/(1-2*ν)


function run()

    # deformation gradient

    F(∇u) = one(∇u) + ∇u
    J(F) = det(F)
    C(F) = (F') * F

    # body force
    B(x) = VectorValue(
        0.0, -0.5, 0.0
    )

    # surfaceTraction
    T(x) = VectorValue(
        -0.1, 0.0, 0.0
    )

    # constitutive relation
    @law function S(∇u)
        Js = J(F(∇u))
        μ * (F(∇u) - (inv(F(∇u)))') + λ * Js * (Js - 1) * (inv(F(∇u)))' 
    end

    function res(u, v)
        S(∇(u)) ⊙ ∇(v) - B ⊙ v
    end

    @law function Lv(∇du, ∇u)
        Js = J(F(∇u))
        Finv = inv(F(∇u))
        μ * (∇du + Finv' ⋅ ∇du' ⋅ Finv') + λ * (2. * Js - 1.) * Js * ((Finv' ⊙ ∇du) * Finv') - λ * Js * (Js - 1) * (Finv' ⋅ ∇du' ⋅ Finv') 
    end

    jac(u, du, v) = Lv(∇(du), ∇(u)) ⊙ ∇(v)

    domain = (0, 1, 0, 1, 0, 1)
    partition = (24, 16, 16)
    tmodel = CartesianDiscreteModel(domain, partition)
    model = simplexify(tmodel)
    labels = get_face_labeling(model)
    add_tag_from_tags!(labels,"rightFace",[2, 4, 6, 8, 14, 16, 18, 20, 26])
    add_tag_from_tags!(labels,"leftFace",[1, 3, 5, 7, 13, 15, 17, 19, 25])
    # add_tag_from_tags!(labels,"topFace",[4, 24])

    V = TestFESpace(
    model=model,valuetype=VectorValue{3,Float64},
    reffe=:Lagrangian,conformity=:H1,order=1,
    dirichlet_tags = ["leftFace", "rightFace"])

    trian = Triangulation(model)
    degree=2
    quad=CellQuadrature(trian, degree)

    neumannTags = ["boundary"]
    btrian = BoundaryTriangulation(model, neumannTags)
    bquad = CellQuadrature(btrian, degree)

    a_Ω = FETerm(res, jac, trian, quad)

    b_Γ(v) =  T ⊙ v
    t_Γ = FESource(b_Γ, btrian, bquad)

    nls = NLSolver(
    show_trace=true,
    method=:newton,
    linesearch=BackTracking()
    )

    g0(x) = VectorValue(0.0,0.0,0.0)
    g1(x) = VectorValue(
        0.0,
        (0.5 + (x[2]-0.5) * cos(π/3) - (x[3] - 0.5) * sin(π/3) - x[2])/2,
        (0.5 + (x[2]-0.5) * sin(π/3) + (x[3] - 0.5) * cos(π/3) - x[3])/2
    )

    U = TrialFESpace(V, [g0, g1])

    #FE problem
    op = FEOperator(U, V, a_Ω, t_Γ)
    solver = FESolver(nls)

    x0 = zeros(Float64, num_free_dofs(V))
    uh = FEFunction(U, x0)
    uh, = solve!(uh, solver, op)

    writevtk(trian,"results_hyper",cellfields=["uh"=>uh,"sigma"=>S(∇(uh))])
end
@btime run()

Output (with v0.13.4)

Iter     f(x) inf-norm    Step 2-norm
------   --------------   --------------
     0     1.078283e-01              NaN
     1     1.707144e-02     8.525316e+01
     2     2.335211e-03     1.941730e+00
     3     1.274437e-04     5.755156e-02
     4     3.548430e-07     1.639837e-05
     5     4.613813e-12     7.484552e-11
  17.212 s (22622611 allocations: 3.12 GiB)

[fverdugo]

@bhaveshshrimali I have also run into this problem, but unfortunately I was not able to create a reproducible MWE.

I have even opened an issue at JuliaLang/julia#37176 but I had to close it since my MWE was not reproducible.

Strangely, if you run (1,) .== 1 at the beginning of your session, the error is fixed. I guess that the error is not related with Gridap, but who knows...

It would be great if you can help to find a reproducible MWE. Then, we can reopen the issue at JuliaLang/julia#37176

[bhaveshshrimali]

Thanks @fverdugo for a prompt response.

I can try some things in a clean REPL tomorrow to see if I can create a reproducible MWE.
Weirdly enough, as you say simply doing (1,).==1 fixes it and I'm now able to run with v0.14.0 absolutely fine without any issues. This is some black magic I have seen today

[fverdugo]

Hope we can find a reproducible MWE and then ask for help to the Julia developers!

[bhaveshshrimali]

Unfortunately, I am not able to reproduce this in a clean REPL. I did a fresh install Julia 1.5.1 on my laptop, and then installed Gridap. Everything works as it should without any problems.
I can close this for now, and it can be reopened if someone stumbles upon this.
@fverdugo ?