I have implemented my own `Type`

of `Real`

(let’s call it `My`

)

Julia’s boast holds true: having defined “all” the operations for `My`

, Julia is “unreasonably” effective at giving me all sorts of matrix algebra operations at zero additional work on my side.

But then I want to compute the `LinearAlgebra.norm`

, and things go sour, although I have implemented all the maths I expect to be needed for an Euclidean norm.

What I do not understand is *what do I need to provide* to get unreasonnable success.

My attempt below is problematic because

- I don’t get a norm
- Error messages have pushed me into implementing conversion to
`Float64`

. In the real world (as opposed to this MnonWE) this conversion would be lossy, and my norm would be wrong
- Error messages have pushed me to implement
`LinearAlgebra.norm`

. Here I have lost the"unreasonnable effectivity" of not having to think what packages `My`

is going to compile together with.

Of course I could write my own `norm`

function, but then I loose any `LinearAlgebra`

functionality that needs a norm…

So here are my attempts.

```
struct My <: Real
x::Float64
end
Base.:(+)(a::My,b::My) = My(a.x+b.x)
Base.:(-)(a::My,b::My) = My(a.x-b.x)
Base.:(*)(a::My,b::My) = My(a.x*b.x)
Base.:(/)(a::My,b::My) = My(a.x/b.x)
Base.:(^)(a::My,b::My) = My(a.x^b.x)
Base.:(^)(a::My,b::Real) = My(a.x^b)
Base.abs(a::My) = My(abs(a.x))
Base.sqrt(a::My) = My(sqrt(a.x))
Base.isnan(a::My) = isnan(a.x)
Base.:(<)(a::My,b::My) = a.x<b.x
Base.AbstractFloat(a::My) = a.x
Base.convert(::Type{Float64},a::My) = a.x
using LinearAlgebra
LinearAlgebra.norm(a::My) = abs(a)
v = My.([1,2,3])
@show norm(v)
```

You need to define promotion rules. Right now, `+(::My, ::My)`

works, but `+(::My, ::AnotherNumberType)`

and `+(::AnotherNumberType, ::My)`

do not. In your example, `norm`

fails because it calls `+(x::Float64, y::My)`

, which has a fallback that promotes its arguments to a common type. To handle this, you need to define what common type that should be. This is done by adding `promote_rule`

methods. One example could be

```
Base.promote_rule(::Type{My}, ::Type{T}) where {T<:Number} = promote_type(Float64, T)
```

That way, the type `My`

behaves like a `Float64`

when promoted.

You don’t need to implement *conversion* to `Float64`

, you only need to implement *promotion* with `Float64`

. e.g.

```
Base.promote_rule(::Type{My}, ::Type{Float64}) = My
```

(The reason for this is that `norm`

needs to promote any type to a floating-point type for accumulation, and for accuracy/overflow reasons asks you to promote to at least double precision.)

You’ll also want:

```
Base.float(x::My) = x
```

(`My`

is already a floating-point type … more generally, if you had `My{T}`

then you want to convert to `My{float(T)})`

.

2 Likes

Yes, I have a MWE now, thanks to your advice. The norm is of type `My`

, so that resolves my concern about a lossy conversion. I had to implement `Base.isnan`

, `Base.inf`

and as you said `Base.float`

, and no need for overloading `LinearAlgebra.norm`

(good!!!)

I had to make another change though:

```
struct My <: Number
```

instead of

```
struct My <: Real
```

Now I wonder how my *real* code is going to like that! - let’s find out! … YES

In practice you’d probably want a `My`

type out of promotion. e.g. suppose `My`

is a quaternion type.

More generally, you probably want something like `My`

to be a parameterized type:

```
struct My{T<:Real} <: Real
x::T
end
Base.:(+)(a::My,b::My) = My(a.x+b.x)
Base.:(-)(a::My,b::My) = My(a.x-b.x)
Base.:(*)(a::My,b::My) = My(a.x*b.x)
Base.:(/)(a::My,b::My) = My(a.x/b.x)
Base.:(^)(a::My,b::My) = My(a.x^b.x)
Base.:(^)(a::My,b::Real) = My(a.x^b)
Base.abs(a::My) = My(abs(a.x))
Base.sqrt(a::My) = My(sqrt(a.x))
Base.isnan(a::My) = isnan(a.x)
Base.isinf(a::My) = isinf(a.x)
Base.isfinite(a::My) = isfinite(a.x)
Base.:(<)(a::My,b::My) = a.x<b.x
Base.promote_rule(::Type{My{T}}, ::Type{S}) where {T,S<:Real} = My{promote_type(Float64, T)}
Base.float(x::My) = My(float(x.x))
Base.convert(::Type{My{T}}, x::Real) where {T} = My(convert(T, x))
Base.convert(::Type{My{T}}, x::My) where {T} = My(convert(T, x.x))
```

and then you get

```
julia> using LinearAlgebra
julia> norm(My.([1,2,3]))
My{Float64}(3.7416573867739413)
```

1 Like

Note also that it’s easier to write this sort of code as:

```
for f in (:abs, :sqrt)
@eval Base.$f(a::My) = My($f(a))
end
for f in (:isnan, :isinf, :isfinite)
@eval Base.$f(a::My) = $f(a)
end
```

1 Like

Yes!

My real code now works, thanks to you, and it has all the features you describe: parametric type, eval and convert.

But interestingly, you guys seem to get it to work with `My <: Real`

while I had to move to `My<:Number`

(well, if it works…)

To get everything to work, inclusive `LinearAlgebra.vecnormalise`

I ended up overloading `Base.isinf, isginite, typemax, typemin, floatmax, floatmin, eps`

and `float`

.