Могут ли хорошие системы типов различать матрицы в разных базисах?

Question

Моя программа (Хартри-Фок / итеративный SCF) имеет две матрицы F и F ', которые на самом деле являются одной и той же матрицей, выраженной в двух разных базисах. Я просто потерял три часа времени отладки, потому что я случайно использовал F 'вместо F. В C ++ средство проверки типов не улавливает ошибки такого рода, потому что обе переменные являются Eigen::Matrix<double, 2, 2> объектами.

Мне было интересно, для Haskell / ML / и т.д. люди, если бы вы писали эту программу, вы бы построили систему типов, где F и F 'имели разные типы? Как это будет выглядеть? Я в основном пытаюсь понять, как я могу перенести некоторые логические ошибки на проверку типов.

Редактировать: основа матрицы похожа на единицу. Вы можете сказать 1 л или сколько угодно галлонов, они оба означают одно и то же. Или, чтобы привести пример вектора, вы можете сказать (0,1) в декартовых координатах или (1, pi / 2) в полярных. Но даже если значение одно и то же, числовые значения различны.

Редактировать: Возможно, единицы были неправильной аналогией. Я не ищу какой-то тип записи, где я могу указать, что первое поле будет литрами, а вторые галлоны, а скорее способ сказать, что эта матрица в целом определяется в терминах какой-то другой матрицы ( базис), где базисом может быть любая матрица одинаковых размеров. Например, конструктор будет выглядеть примерно как mkMatrix [[1, 2], [3, 4]] [[5, 6], [7, 8]], а затем добавление этого объекта в другую матрицу будет проверять тип, только если у обоих объектов будет та же матрица, что и у их вторых параметров. Имеет ли это смысл?

Редактировать: определение Википедия , рабочие примеры

Answer 1

Это вполне возможно в Хаскеле.

Статически проверенные размеры

Haskell имеет массивы с статически проверенными размерами , где измерения можно манипулировать и проверять статически, предотвращая индексацию в неправильном измерении. Некоторые примеры:

Это будет работать только на двумерных массивах :

multiplyMM :: Array DIM2 Double -> Array DIM2 Double -> Array DIM2 Double

Пример из repa должен дать вам смысл. Здесь, для диагонали требуется двумерный массив, возвращается одномерный массив того же типа.

diagonal :: Array DIM2 e -> Array DIM1 e

или, из Руководство по ремонту Matt sottile , статически проверенные размеры в трехмерном матричном преобразовании:

f :: Array DIM3 Double -> Array DIM2 Double
f u =
  let slabX = (Z:.All:.All:.(0::Int))
      slabY = (Z:.All:.All:.(1::Int))
      u' = (slice u slabX) * (slice u slabX) +
           (slice u slabY) * (slice u slabY)
  in
    R.map sqrt u'

Статически проверенные единицы

Другой пример вне матричного программирования: статически проверено единицы измерения , что приводит к ошибке типа, чтобы запутать, например. футы и метры, не делая преобразования.

 Prelude> 3 *~ foot + 1 *~ metre
 1.9144 m

или для всего набора единиц СИ и количеств .

например. Нельзя добавлять вещи другого измерения, например, объемы и длины:

> 1 *~ centi litre + 2 *~ inch 
Error:
Expected type: Unit DVolume a1
  Actual type: Unit DLength a0

Итак, следуя типам измерений массива в стиле repa , я бы предложил добавить параметр фантомного типа Base к типу массива и использовать его для различения баз. В Хаскеле индекс Dim Аргумент типа дает ранг массива (то есть его форму), и вы можете сделать то же самое.

Или, если под базой вы подразумеваете какое-то измерение в единицах, используя типы измерений.

Итак, да, сейчас это почти обычная техника в Haskell, и есть несколько примеров проектирования с такими типами, чтобы помочь вам начать.

Answer 2

Это очень хороший вопрос. Я не думаю, что вы можете закодировать понятие базиса в большинстве систем типов, потому что по сути все, что делает средство проверки типов, должно иметь возможность завершаться, и делать суждения о том, равны ли два действительных вектора, слишком сложно. Вы можете иметь (2 v_1) + (2 v_2) или 2 (v_1 + v_2), например. Есть некоторые языки, которые используют зависимые типы [ wikipedia ], но они относительно академичны.

Я думаю, что большая часть вашей отладочной боли будет уменьшена, если вы просто закодируете основы, в которых ваша матрица работает вместе с матрицей. Например,

newtype Matrix = Matrix { transform :: [[Double]], 
    srcbasis :: [Double], dstbasis :: [Double] }

и затем, когда вы M из базы a до b с N , проверьте, что N от b до c , и вернуть матрицу с базой a до c .

ПРИМЕЧАНИЕ - кажется, что большинство людей здесь имеют программирование вместо математического фона, поэтому я приведу краткое объяснение здесь. Матрицы являются кодировками линейных преобразований между векторными пространствами. Например, если вы кодируете поворот на 45 градусов в R ^ 2 (2-мерные числа), то стандартный способ кодирования этого в матрице состоит в том, что стандартный базис вектор e_1, записанный «[1, 0]», отправляется на комбинацию e_1 и e_2, а именно [1 / sqrt (2), 1 / sqrt (2)]. Дело в том, что вы можете закодировать один и тот же поворот, сказав, куда идут разные векторы, например, вы можете сказать, куда вы отправляете [1,1] и [1, -1] вместо e_1 = [1,0] и e_2 = [0,1], и это будет иметь другое матричное представление.

Редактировать 1

Если у вас есть конечный набор баз, с которыми вы работаете, вы можете сделать это ...

{-# LANGUAGE EmptyDataDecls #-}
data BasisA
data BasisB
data BasisC

newtype Matrix a b = Matrix { coefficients :: [[Double]] }

multiply :: Matrix a b -> Matrix b c -> Matrix a c
multiply (Matrix a_coeff) (Matrix b_coeff) = (Matrix multiplied) :: Matrix a c
    where multiplied = undefined -- your algorithm here

Затем, в ghci (интерактивный интерпретатор Haskell),

*Matrix> let m = Matrix [[1, 2], [3, 4]] :: Matrix BasisA BasisB
*Matrix> m `multiply` m

<interactive>:1:13:
    Couldn't match expected type `BasisB'
        against inferred type `BasisA'
*Matrix> let m2 = Matrix [[1, 2], [3, 4]] :: Matrix BasisB BasisC
*Matrix> m `multiply` m2
-- works after you finish defining show and the multiplication algorithm

Answer 3

Хотя я понимаю, что это не решает строго (уточненный) вопрос - мои извинения - оно кажется актуальным, по крайней мере, в отношении популярного ответа Дона Стюарта ...

Я являюсь автором библиотеки Haskell dimension , на которую ссылался Дон и приводил примеры. Я также писал - в некоторой степени под радаром - экспериментальную элементарную библиотеку линейной алгебры, основанную на размерной . Эта библиотека линейной алгебры статически отслеживает размеры векторов и матриц, а также физические размеры («единицы») их элементов на основе каждого элемента.

Эта последняя точка - отслеживание физических измерений для каждого элемента - является довольно сложной и, возможно, избыточной для большинства применений, и можно даже утверждать, что не имеет большого математического смысла иметь количества различных физических измерений в качестве элементов в любом данном векторе / матрица. Однако некоторые приложения линейной алгебры, представляющие интерес для меня, такие как фильтрация Калмана и взвешенная оценка наименьших квадратов, обычно используют гетерогенные векторы состояний и ковариационные матрицы.

Используя фильтр Калмана в качестве примера, рассмотрим вектор состояния x = [d, v] , который имеет физические размеры [L, LT ^ -1] . Следующий (будущий) вектор состояния прогнозируется путем умножения на матрицу перехода состояний F , то есть: x '= F x_ . Очевидно, что это уравнение имеет смысл F не может быть произвольным, но должно иметь размер и физические размеры [[1, T], [T ^ -1, 1]] . Функция predict_x', представленная ниже, статически обеспечивает выполнение этого отношения:

predict_x' :: (Num a, MatrixVector f x x) => Mat f a -> Vec x a -> Vec x a
predict_x' f x_ = f |*< x_

(неприглядный оператор |*< обозначает умножение матрицы слева на вектор справа.)

В более общем случае для априорного вектора состояния x_ произвольного размера и с элементами произвольных физических измерений передача матрицы перехода состояний f с "несовместимым" размером и / или физическими размерами в predict_x' вызовет ошибка времени компиляции.

Answer 4

В F # (который первоначально произошел от OCaml), вы можете использовать единицы измерения .Эндрю Кеннед, который разработал эту функцию (а также создал очень интересную теорию), написал серию статей , демонстрирующих ее .

. Это вполне может быть использовано в вашем сценарии -хотя я не до конца понимаю вопрос.Например, вы можете объявить два типа единиц измерения следующим образом:

[<Measure>] type litre
[<Measure>] type gallon

Добавление литров и галлонов дает ошибку времени компиляции:

1.0<litre> + 1.0<gallon> // Error!

F # не вставляет автоматически преобразование между различнымиединиц, но вы можете написать функцию преобразования:

let toLitres gal = gal * 3.78541178<litre/gallon>
1.0<litre> + (toLitres 1.0<gallon>)

Прекрасная вещь в единицах измерения в F # состоит в том, что они автоматически выводятся, а функции являются общими.Если вы умножите 1.0<gallon> * 1.0<gallon>, результат будет 1.0<gallon^2>.

Люди использовали эту функцию для разных целей - от преобразования виртуальных счетчиков до пикселей экрана (при моделировании солнечной системы) до преобразования валют (долларов вфинансовые системы).Несмотря на то, что я не эксперт, вполне вероятно, что вы можете каким-то образом использовать его и для проблемной области.

Answer 5

Если это выражено в другой базе, вы можете просто добавить параметр шаблона, который будет действовать как база. Это будет дифференцировать эти типы. float - это float - это float - если вы не хотите, чтобы два значения float были одинаковыми, если они на самом деле имеют одно и то же значение, вам нужно сообщить об этом системе типов.