CliMA / ClimaCore.jl

CliMA model dycore
https://clima.github.io/ClimaCore.jl/dev
Apache License 2.0
87 stars 8 forks source link

Add debug compare #711

Open charleskawczynski opened 2 years ago

charleskawczynski commented 2 years ago

I think we should add support for a highly generic debug_compare function for users to compare two Fields / FieldVectors. I often find myself comparing two fields / field vectors, and I often want to know: 1) do the two Fields/FieldVectors have the same field names and 2) if they do, what is the difference in their values.

Describe the solution you'd like

We have something like this in TurbulenceConvection.jl which looks something like this:

"""
    debug_compare(a::Field, a::Field)
    debug_compare(a::FieldVector, b::FieldVector)

Recursively compare two identically structured
`Field`s, or `FieldVector`s, with potentially different
data. If `!(maximum(abs.(parent(a) .- parent(b))) == 0.0)`
for any single field, then `debug_compare` will print out the fields
and `err`. This can be helpful for debugging where and why
two `Field`/`FieldVector`s are different.
"""
function debug_compare(a::FVA, b::FVB, pn0 = "") where {FVA <: ClimaCore.Fields.FieldVector, FVB <: ClimaCore.Fields.FieldVector}
    for pn in propertynames(a)
        pa = getproperty(a, pn)
        pb = getproperty(b, pn)
        debug_compare(pa, pb, "$pn0.$pn")
    end
end

function debug_compare(a::F, b::G, pn0 = "") where {F <: ClimaCore.Fields.Field, G <: ClimaCore.Fields.Field}
    if isempty(propertynames(a))
        if any(size(parent(a)) .≠ size(parent(b)))
            @info "Field `$pn0` has non-comparable size: size(a)=$(size(parent(a))),  size(b)=$(size(parent(b)))"
        else
            err = abs.(parent(a) .- parent(b))
            if !(maximum(err) == 0.0)
                println("--- Comparing field $pn0")
                @show a
                @show b
                @show err
            end
        end
    else
        for pn in propertynames(a)
            pa = getproperty(a, pn)
            pb = getproperty(b, pn)
            debug_compare(pa, pb, "$pn0.$pn")
        end
    end
end

Can we add this? cc @simonbyrne, @valeriabarra, @dennisYatunin

I was thinking that we could add verbosity options to print more or less info.

valeriabarra commented 2 years ago

I can see this could be a useful feature. Maybe the name could just be compare without the debug bit? In case people can use this for any other type of investigation, Idk...