Решение той задачи, что вы привели, в аналитическом виде скорее всего будет очень-очень сложным. Напряженность снаружи соленоида в любой точке - это еще полдела (хотя и это чрезвычайно сложно, и я не уверен, что это вообще кто-то решал), настоящие проблемы вылезут, когда движущийся магнит будет создавать в катушке противоЭДС, которым в двигателе пренебрегать никак нельзя. Боюсь представить, как это все будет выглядеть. И даже если каким-то невероятным образом это удастся решить, к пониманию работы BLDC вряд ли приведет из-за сложности итоговых уравнений.
Я бы рассмотрел двигатель как магнитную цепь, вывел бы из нее значения напряженности магнитного поля статора и магнитного потока в роторе (зависящие от его положения и скорости), откуда легко считается момент сил, действующий на ротор. Решение будет очень приближенным, но его вполне можно будет использовать, особенно если померить реальные характеристики расчетного двигателя и внести на глаз в аналитическое решение какие-нибудь поправочные коэффициенты или нелинейные множители. Приближенность решения обуславливается допущениями, принимаемыми теорией магнитных цепей. Если вы с ней не знакомы (или подзабыли), я попробую вкратце про нее написать под оффтопом (хоть и тоже последний раз пользовался ей лет 15 назад).
(Оффтоп)
В магнитной цепи рассматривается "течение" магнитного потока по магнитопроводам и воздушным зазорам между ними. Проводятся аналогии с течением электрического тока в проводниках, строится электрическая схема замещения магнитной системы, затем она рассчитывается как электрическая схема, и таким образом находятся магнитные потоки на всех участках. При этом считается, что (все формулы в системе СИ):
1. Магнитопровод не намагничивается, напряженность магнитного поля в нем линейно зависит от магнитной индукции (

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

, где B - магнитная индукция, S - поперечное сечение магнитопровода), который "течет" по магнитопроводу.
3. В воздушных зазорах между магнитопроводами на пути "течения" потока поле также однородно и не выходит за пределы, ограниченные минимальным поперечным сечением магнитопроводов по соседству (понятно будет на картинке ниже) - при достаточно малых зазорах это допущение не приводит к особым неточностям.
Определения:
1. Магнитодвижущая сила (МДС) - аналог ЭДС из электрических цепей - величина, характеризующая намотанную вокруг магнитопровода катушку с током, создающую магнитный поток (как ЭДС характеризует "батарейку", создающую ток). Равна

, где I - ток в катушке, n - число витков. Направление МДС соответствует направлению создаваемого потока и определяется по правилу буравчика.
2. Магнитное сопротивление - аналог электрического сопротивления - величина, характеризующая "препятствие" магнитопровода к "течению" магнитного потока. Равна

, где l - длина магнитопровода, S - площадь его поперечного сечения,

- магнитная проницаемость материала магнитопровода,

- магнитная постоянная. Для воздушных зазоров

. Сопротивления последовательно соединенных магнитопроводов складываются.
Тут еще важно понимать, что в отличие от электрических цепей "протекание" постоянного потока по магнитного сопротивлению не вызывает никаких энергетических потерь, то есть аналогия чисто математическая.
3. Магнитное напряжение - аналог электрического напряжения - "разность потенциалов" магнитного поля на какой-то длине магнитопровода (или воздушного зазора). Отдельной буквой не обозначается, равна

, где H - напряженность магнитного поля, l - длина участка магнитопровода (или зазора), на котором определяется напряжение.
Первый закон Кирхгофа - сумма потоков, втекающих в узел, равна сумме потоков, вытекающих из узла.
Второй закон Кирхгофа - сумма МДС в контуре равна сумме магнитных напряжений в нем.
Закон Ома для участка цепи -

В принципе, в интернете много информации по магнитным цепям. Алгоритм их построения и расчета такой:
1. Разбиваете весь магнитопровод на участки, для которых можно выделить длину и среднюю площадь. Выделяете воздушные зазоры.
2. Считаете для каждого участка и зазора магнитные сопротивления.
3. Считаете все МДС.
4. Добавляете МДС и сопротивления на схему, соединяя их так, как соответствующие объекты стыкуются в реальности.
5. Рассчитываете магнитные потоки и напряжения как для электрической схемы, используя закон Ома и законы Кирхгофа.
Постоянный магнит можно рассматривать как "источник магнитного потока" по аналогии с источником тока в электрической цепи.
В вашем случае магнитная цепь будет с переменным воздушным зазором, поскольку ротор вращается относительно статора. Длина зазора будет постоянной, а меняться будет площадь. Соответственно, будет меняться и магнитное сопротивление - его надо выразить в виде функции от угла положения ротора.
Из-за переменного сопротивления будет меняться и магнитный поток. По закону электромагнитной индукции в катушках будет возникать ЭДС, препятствующая этому изменению. Это тоже обязательно надо учесть в системе уравнений, иначе магнитный поток получится завышенным.
В итоге из расчета этой схемы нужно вывести выражение для магнитного потока через ротор и напряженность поля статора. Тогда момент магнитных сил, действующий на ротор, будет равен

, где d - диаметр ротора,

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

, где

- момент нагрузки, J - момент инерции ротора,

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