/* Соглашения о нумерации узлов. nd - количество диодов. 1..nd+1 - последовательно пронумерованные узлы диодной цепи. 1 - первый узел диодной цепи, является общим проводом в случае чётного коэффициента; подключается к источнику в случае нечётного. nd+1 - последний, выходной узел диодной цепи, к нему подключается нагрузка. Номер узла, соответствующего общему проводу, должен обязательно быть равен 1 (чётный коэффициент), либо иметь наибольший номер среди всех узлов. Первый вывод источника питания подключается ко входу умножителя: к узлу 1, если умножитель имеет нечётный коэффициент или к некоторому конденсатору, если коэффициент чётный. Второй вывод источника питания подключается всегда к общему проводу (так можно определить номер узла, соответствующего общему проводу). Если общий провод имеет номер 1, то максимальный номер должен быть у узла- входа умножителя (к которому подключён источник). В результате принятых соглашений, питающий источник оказывается всегда включённым между узлами 1 и узлом с наибольшим номером. */ /* Соглашения о структуре матрицы соединений circuit. Если схема умножителя содержит nc конденсаторов, матрица будет иметь размер (nc+1)*2. Первые строки 1..nc описывают подключение соответствующих конденсаторов: circuit[i, 1] - номер узла, к которому подключён первый вывод конденсатора Ci; circuit[i, 2] - номер узла, к которому подключён второй вывод. Последняя строка описывает подключение источника питания. Диоды и нагрузка не указываются в матрице, способ их подключения является предопределённым. */ /* Присвоить элементам строки i матрицы A значения из списка (матрицы) r */ set_row(A, i, r):=block( [], if not listp(r) then r:list_matrix_entries(r), for j:1 thru length(r) do A[i, j]:r[j], A)$ /* Присвоить элементам столбца j матрицы A значения из списка (матрицы) c */ set_col(A, j, c):=block( [cl], cl:if not listp(c) then list_matrix_entries(c) else c, for i:1 thru length(cl) do A[i, j]:cl[i], A)$ /* Функция возвращает матрицу соединений для умножителя на nd типа A. */ create_circuitA(nd):=block( [n, m, circuit], n:nd+2, m:nd+1, circuit:zeromatrix(m, 2), for i:1 thru nd do ( circuit[i, 1]:i+1, circuit[i, 2]:i-1 ), circuit[1, 2]:n, circuit[m, 1]:if mod(nd, 2)=0 then n else 1, circuit[m, 2]:if mod(nd, 2)=0 then 1 else n, circuit)$ /* Функция возвращает матрицу соединений для умножителя на nd типа B. */ create_circuitB(nd):=block( [n, m, circuit], if nd<=2 then return(create_circuitA(nd)), n:nd+2, m:nd+1, circuit:zeromatrix(m, 2), for i:1 thru m do ( circuit[i, 1]:i+1, circuit[i, 2]:i-3 ), circuit[1, 2]:n, circuit[2, 2]:1, circuit[3, 2]:n, circuit[m, 2]:1, if mod(nd, 2)=1 then ( circuit[m, 1]:1, circuit[m, 2]:n ), circuit)$ /* Функция возвращает матрицу соединений для умножителя на nd типа Z. */ create_circuitZ(nd):=block( [n, m, circuit], n:nd+2, m:nd+1, circuit:zeromatrix(m, 2), for i:1 thru m do ( circuit[i, 1]:i+1, circuit[i, 2]:if mod(i, 2)=1 then n else 1 ), circuit[m, 2]:1, if mod(nd, 2)=1 then ( circuit[m, 1]:1, circuit[m, 2]:n ), circuit)$ /* По матрице соединений сформировать узловую матрицу, считая диоды разомкнутыми ветвями (т.е. игнорируя их). Автоматически добавляется ветвь - нагрузка. */ get_node_matrix(circuit, nd):=block( [nn, nbr, gnd_node, A], /* TODO: можно добавить проверку валидности circuit */ /* Ищем узел с наибольшим номером, считаем этот номер равным общему количеству узлов. */ nn:lmax(list_matrix_entries(circuit)), /* Количество ветвей равно количеству строк в circuit + 1 ветвь - нагрузка */ nbr:length(circuit)+1, A:zeromatrix(nn, nbr), for i:1 thru nbr-1 do ( A[circuit[i, 1], i]:1, A[circuit[i, 2], i]:-1), gnd_node:circuit[nbr-1, 2], A[nd+1, nbr]:1, A[gnd_node, nbr]:-1, A)$ /* Возвращает узловую матрицу, полученную из узловой матрицы A исключением висящих ветвей. */ exclude_orphan_branches(A):=block( [A:copymatrix(A), sz, r, k, j, n:1], sz:matrix_size(A), /* Выполняем в цикле, т.к. исключение одних висящих ветвей может сделать висящими другие ветви. */ while n#0 do ( n:0, for i:1 thru sz[1] do ( r:abs(row(A, i)), if sum(r[1, k], k, 1, sz[2])=1 then ( j:locate_matrix_entry(r, 1, 1, 1, sz[2], identity, 'max)[2], set_col(A, j, zeromatrix(sz[1], 1)), n:n+1))), A)$ /* Возвращает узловую матрицу дерева, построенного для схемы, которая задана узловой матрицей A. */ get_tree_matrix(A):=block( [T, w, sz, c, z, ntr, Test, s], T:zerofor(A), sz:matrix_size(A), w:zeromatrix(sz[1], 1), ntr:{}, for j:1 thru sz[2] do ( c:abs(col(A, j)), z:w*c, s:sum(z[i, 1], i, 1, length(z)), if s<2 then ( w:w+c-z, set_col(T, j, col(A, j)) ) elseif s=2 then ( Test:copymatrix(T), set_col(Test, j, col(A, j)), if zeromatrixp(exclude_orphan_branches(Test)) then ( set_col(T, j, col(A, j)) ) else ntr:adjoin(j, ntr) ) ), [T, ntr])$ /* Возвращает строку матрицы контуров, соответствующую одному контуру, найденному по узловой матрице A. Возвращаемый вектор-строка имеет тип список. */ get_closed_loop(A):=block( [A:copymatrix(A), n, m, sz, result, nl, ie, i, j, r], [n, m]:matrix_size(A), nl:makelist(0, m), result:copylist(nl), A:exclude_orphan_branches(A), [i, j]:locate_matrix_entry(A, 1, 1, n, m, lambda([x], x#0), 'bool), if i=false then return(nl), result[ j]:A[i, j], A[i, j]:0, [ie, j]:locate_matrix_entry(A, 1, j, n, j, lambda([x], x#0), 'bool), if ie=false then return(nl), A[ie, j]:0, r:while i#ie do ( [i, j]:locate_matrix_entry(A, i, 1, i, m, lambda([x], x#0), 'bool), if j=false then return(false), result[j]:-A[i, j], A[i, j]:0, [i, j]:locate_matrix_entry(A, 1, j, n, j, lambda([x], x#0), 'bool), if j=false then return(false), A[i, j]:0 ), if r=false then return(nl), result)$ /* Возвращает матрицу контуров для схемы, заданной узловой матрицей A. */ get_independent_loops(A):=block( [T, ntr, lst, result, r, tmp], [T, ntr]:get_tree_matrix(A), lst:listify(ntr), result:matrix(), for j:1 thru length(lst) do ( tmp:copymatrix(T), set_col(tmp, lst[j], col(A, lst[j])), r:get_closed_loop(tmp), if not zeromatrixp(r) then result:addrow(result, r) ), result)$ /* Возвращает узловую матрицу, получаемую из A при условии замыкания нечётных (k=1) или чётных (k=2) диодов. */ switch_diodes(A, nd, k):=block( [result, n, m, n2, nj, r1, r2, di, i0], if k#1 and k#2 then error("Недопустимое значение аргумента k."), [n, m]:matrix_size(A), nj:if mod(nd, 2)=1 then (nd+3-2*k)/2 else nd/2, n2:n-nj, result:zeromatrix(n2, m), di:if k=1 then 0 else 1, if di>0 then set_row(result, 1, row(A, 1)), for i:1 thru nj do ( r1:row(A, 2*i-1+di), r2:row(A, 2*i+di), if not zeromatrixp(r1*r2) then error("Недопустимая узловая матрица: конденсатор включён между соседними узлами диодной цепи."), set_row(result, i+di, r1+r2) ), i0:2*nj+di+1, for i:i0 thru n do set_row(result, i-nj, row(A, i)), result)$ /* Возвращает уравнения переключения - уравнения цепи (в виде списка уравнений), составленные по узловой матрице A при условии замыкания нечётных (k=1) или чётных диодов (k=2). nd - количество диодов в умножителе. x1, x2 - список переменных в начальный и конечный моменты времени. xc - список с ёмкостями конденсаторов. */ get_switch_equations(A, nd, k, x1, x2, xc):=block( [As, An, Y, X, eq, n, m, Q], As:switch_diodes(A, nd, k), /* Узловые уравнения */ An:submatrix(1, length(As), As), [n, m]:matrix_size(An), X:if m=0 then [] else append(makelist(xc[i]*(x2[i]-x1[i]), i, 1, m-2), [0, 0]), Y:An.X, if not matrixp(Y) then Y:matrix([Y]), eq:makelist(Y[i, 1]=0, i, 1, n), /* Контурные уравнения */ Q:get_independent_loops(As), [n, m]:matrix_size(Q), if n=0 then return(eq), X:makelist(x2[i], i, 1, m), Y:Q.X, if not matrixp(Y) then Y:matrix([Y]), n:length(Q), eq:append(eq, makelist(Y[i, 1]=0, i, 1, n)), eq)$ /* Функция составляет систему уравнений для заданного умножителя и решает её. Возвращает список из трёх элементов: выходное напряжение, размах пульсаций и список из всех найденных решений системы уравнений (напряжения на конденсаторах, мгновенное напряжение на источнике и на выходе в характерных точках). */ solve_multiplier(circuit, nd, xc, a, F, IL):=block( [u, ux, u1, uy, u2, A, T, nbr, n, m, out_node, gnd_node, power_node, M1, M2, X1, X2, Y, eq, g, k, xlst, r, result, uout, hout], if not listp(xc) or length(xc)#length(circuit)-1 then error("Аргумент xc должен быть списком, количество элементов которого равно количеству конденсаторов в схеме."), A:get_node_matrix(circuit, nd), [n, m]:matrix_size(A), [T, nbr]:get_tree_matrix(A), out_node:nd+1, gnd_node:circuit[length(circuit), 2], power_node:circuit[length(circuit), 1], eq:[], /* Уравнения разряда при закрытых диодах, 2 однотипных системы. */ M1:submatrix(gnd_node, power_node, A), X1:makelist(0, m), X2:makelist(0, m), for i:1 thru m-2 do ( X1[i]:xc[i]*(ux[i]-u[i]), X2[i]:xc[i]*(uy[i]-u1[i]) ), X1[m]:IL/(2*F), X2[m]:IL/(2*F), Y:M1.X1, if not matrixp(Y) then Y:matrix([Y]), /* Получаем уравнений на 2 меньше, чем узлов */ eq:append(eq, makelist(Y[i, 1]=0, i, 1, n-2)), Y:M1.X2, if not matrixp(Y) then Y:matrix([Y]), eq:append(eq, makelist(Y[i, 1]=0, i, 1, n-2)), /* Если есть замкнутые контуры, добавляем уравнения для них */ M2:get_independent_loops(A), g:length(M2), if g>0 then ( X1:makelist(ux[i], i, 1, m), X2:makelist(uy[i], i, 1, m), Y:M2.X1, if not matrixp(Y) then Y:matrix([Y]), eq:append(eq, makelist(Y[i, 1]=0, i, 1, g)), Y:M2.X2, if not matrixp(Y) then Y:matrix([Y]), eq:append(eq, makelist(Y[i, 1]=0, i, 1, g)) ), /* Уравнения первого перезаряда: ux->u1 */ k:if mod(nd, 2)=0 then 1 else 2, X1:makelist(ux[i], i, 1, m), X2:makelist(u1[i], i, 1, m), eq:append(eq, get_switch_equations(A, nd, k, X1, X2, xc)), /* Уравнения второго перезаряда: uy->u2 */ k:if mod(nd, 2)=0 then 2 else 1, X1:makelist(uy[i], i, 1, m), X2:makelist(u2[i], i, 1, m), eq:append(eq, get_switch_equations(A, nd, k, X1, X2, xc)), /* Уравнения стационарного режима: u[i]=u2[i] */ eq:append(eq, makelist(u[i]=u2[i], i, 1, m)), /* Уравнения питающего умножитель источника переменного напряжения. */ eq:append(eq, [u[m-1]=a, ux[m-1]=-a, u1[m-1]=-a, uy[m-1]=a]), /* Решаем полученную систему относительно следующих переменных: */ xlst:append( makelist(u[i], i, 1, m), makelist(ux[i], i, 1, m), makelist(u1[i], i, 1, m), makelist(uy[i], i, 1, m), makelist(u2[i], i, 1, m)), r:solve(eq, xlst), if length(r)=0 then error("Ошибка: не найдено решение системы."), /* Преобразуем результат решения в список значений переменных в порядке их перечисления в xlst, значения u2[i] не включаем в список */ result:makelist(args(r[1][i])[2], i, 1, 4*m), uout:result[m], hout:uout-result[4*m], [expand(uout), expand(hout), expand(result)] )$ /* Тестирование функции solve_multiplier на заведомо известных решениях */ test_solve_multiplier(iterations):=block( [tp:["A", "B", "Z"], type, create_circuit, err:1/100, d, f, fh, fa, fha, fb, fhb, fz, fhz, uout, hout, r, th_uout, th_hout, cir, result:"OK"], d:IL/(2*F*C), fa(m):=m*a-d*( if mod(m, 2)=0 then m^3/6+m^2/8+m/12 else m^3/6+m^2/8-m/6-1/8), fha(m):=d*if mod(m, 2)=0 then (m^2/4+m/2) else (m^2/4+m+3/4), fb(m):=m*a-d*(m^3/24+5*m^2/32+ if mod(m, 4)=0 then 5*m/24 elseif mod(m, 4)=1 then m/48-7/32 elseif mod(m, 4)=2 then m/3+3/8 else /*mod(m, 4)=3*/ 19*m/48+9/32), fhb(m):=d*(m^2/16+ if mod(m, 4)=0 then m/4 elseif mod(m, 4)=1 then 5*m/8+21/16 elseif mod(m, 4)=2 then m/2+3/4 else 3*m/8+5/16), fz(m):=m*a-d*(2*m-2), fhz(m):=2*d, f(type, m):=if type="A" then fa(m) elseif type="B" then fb(m) elseif type="Z" then fz(m) else error(type, "is an unknown type of circuits."), fh(type, m):=if type="A" then fha(m) elseif type="B" then fhb(m) elseif type="Z" then fhz(m) else error(type, "is an unknown type of circuits."), create_circuit(type, m):=if type="A" then create_circuitA(m) elseif type="B" then create_circuitB(m) elseif type="Z" then create_circuitZ(m) else error(type, "is an unknown type of circuits."), prederror:true, display("Type A, x12, C1=4C"), th_uout:12*a-253*d, th_hout:42*d, display(th_uout, th_hout), xc:makelist(C, 12), xc[1]:4*C, [uout, hout, r]:solve_multiplier(create_circuitA(12), 12, xc, a, F, IL), if not equal(th_uout, uout) or not equal(th_hout, hout) then error("Обнаружено расхождение."), display(result), display("Type Design1, x4"), th_uout:4*a-17/3*d, th_hout:2*d, display(th_uout, th_hout), xc:makelist(C, 6), cir:matrix([2, 6], [3, 1], [4, 2], [5, 3], [4, 6], [5, 1], [6, 1]), [uout, hout, r]:solve_multiplier(cir, 4, xc, a, F, IL), if not equal(th_uout, uout) or not equal(th_hout, hout) then error("Обнаружено расхождение."), display(result), for nd:1 thru iterations do ( for j:1 thru length(tp) do ( type:tp[j], th_uout:f(type, nd), th_hout:fh(type, nd), xc:makelist(C, nd), display(type, nd, th_uout, th_hout), [uout, hout, r]:solve_multiplier(create_circuit(type, nd), nd, xc, a, F, IL), display(uout, hout), if not equal(th_uout, uout) or not equal(th_hout, hout) then error("Обнаружено расхождение."), display(result) ) ), true)$