Multimaterial Bug (tags)

[zjwegert]

Hi, I've been working on implementing a 2D multi-material plane stress problem in Gridap using Voigt notation. The given stiffness matrix looks good but when adding tags to specify different materials I'm getting an error. Below is the @law functions for the bilinear form a and the definition for a.

@law function σ(ε,tag)
    if tag == -1 # example
        return zeros(3,3)*vector_ε(ε)
    else
        D = E/(1-ν^2)*Matrix([1 ν 0
                              ν 1 0
                              0 0 (1-ν)/2]);
        return D*vector_ε(ε)
    end
end
@law function vector_ε(ε)
    return [ε[1,1];ε[2,2];2*ε[1,2]]
end

a(u,v) = vector_ε(ε(v)) ⋅ σ(ε(u),tags)

This code should just treat the whole multi-material problem as a single material. However, it produces the following output at the line AffineFEOperator(U, V0, t_Ω):

MethodError: no method matching zero(::Type{Array{Float64,1}})
Closest candidates are:
  zero(!Matched::Type{Dates.DateTime}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\Dates\src\types.jl:404
  zero(!Matched::Type{Pkg.Resolve.VersionWeight}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\Pkg\src\Resolve\versionweights.jl:15
  zero(!Matched::Type{Dates.Time}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\Dates\src\types.jl:406
  ...
in include_string at base\loading.jl:1088
in top-level scope at top_levelset.jl:421
in AffineFEOperator at Gridap\Fw0X0\src\FESpaces\AffineFEOperators.jl:44
in AffineFEOperator at Gridap\Fw0X0\src\FESpaces\AffineFEOperators.jl:31
in collect_cell_matrix_and_vector at Gridap\Fw0X0\src\FESpaces\FETerms.jl:236
in get_cell_matrix_and_vector at Gridap\Fw0X0\src\FESpaces\FETerms.jl:107
in get_cell_matrix at Gridap\Fw0X0\src\FESpaces\FETerms.jl:432
in integrate at Gridap\Fw0X0\src\FESpaces\CellBases.jl:368
in integrate at Gridap\Fw0X0\src\Fields\Integrate.jl:17
in evaluate_field_array at Gridap\Fw0X0\src\Fields\FieldArrays.jl:79 
in kernel_evaluate at Gridap\Fw0X0\src\Fields\FieldArrays.jl:274
in evaluate_field_arrays at Gridap\Fw0X0\src\Fields\FieldArrays.jl:47 
in _evaluate_field_arrays at Gridap\Fw0X0\src\Fields\FieldArrays.jl:52 
in evaluate_field_arrays at Gridap\Fw0X0\src\Fields\FieldArrays.jl:47 
in _evaluate_field_arrays at Gridap\Fw0X0\src\Fields\FieldArrays.jl:57 
in evaluate_field_array at Gridap\Fw0X0\src\Fields\FieldArrays.jl:79 
in kernel_evaluate at Gridap\Fw0X0\src\Fields\FieldArrays.jl:275
in apply at Gridap\Fw0X0\src\Arrays\Apply.jl:61
in apply at Gridap\Fw0X0\src\Arrays\Apply.jl:105 
in Gridap.Arrays.AppliedArray at Gridap\Fw0X0\src\Arrays\Apply.jl:176
in kernel_testitem at Gridap\Fw0X0\src\Arrays\Kernels.jl:179
in kernel_cache at Gridap\Fw0X0\src\Fields\FieldOperations.jl:131
in zeros at base\array.jl:526

No error is produced in the single material code with tags removed:

@law function σ(ε)
    tag = 0
    if tag == -1 # example
        return zeros(3,3)*vector_ε(ε)
    else
        D = E/(1-ν^2)*Matrix([1 ν 0
                              ν 1 0
                              0 0 (1-ν)/2]);
        return D*vector_ε(ε)
    end
end
@law function vector_ε(ε)
    return [ε[1,1];ε[2,2];2*ε[1,2]]
end

a(u,v) = vector_ε(ε(v)) ⋅ σ(ε(u))

This looks like a bug.

[fverdugo]

You cannot use plain Vector objects to represent stresses/strains and Matrix objects for the constitutive tensor in Voigt notation. Use VectorValue and TensorValue instead.

On top of that, I do not believe that you really need to use Voig notation in practice. We already have symmetric tensors that represent stress/strain efficiently.