| f g det2x2 partialDerivative jacobian aN bN newtonMethod x0 y0 h epsilon |
" Функции для вычисления определителя матрицы 2х2 и частной производной "
f := [:x :y | x * x + y * y - 4].
g := [:x :y | x - y].
det2x2 := [:a :b :c :d | a * d - b * c].
partialDerivative := [:func :x :y :h :derivativeType |
derivativeType = 'x'
ifTrue: [(func value: x + h value: y) - (func value: x value: y) / h]
ifFalse: [derivativeType = 'y'
ifTrue: [(func value: x value: y + h) - (func value: x value: y) / h]
ifFalse: [Error signal: 'Invalid derivative type']]].
jacobian := [:func1 :func2 :x :y :h |
| dx1 dy1 dx2 dy2 |
dx1 := partialDerivative valueWithArguments: {func1. x. y. h. 'x'}.
dy1 := partialDerivative valueWithArguments: {func1. x. y. h. 'y'}.
dx2 := partialDerivative valueWithArguments: {func2. x. y. h. 'x'}.
dy2 := partialDerivative valueWithArguments: {func2. x. y. h. 'y'}.
^ det2x2 valueWithArguments: {dx1. dy1. dx2. dy2}.
].
aN := [:func1 :func2 :x :y :h |
| fx fy dfy1 dfy2 |
fx := func1 value: x value: y.
fy := func2 value: x value: y.
dfy1 := partialDerivative valueWithArguments: {func1. x. y. h. 'y'}.
dfy2 := partialDerivative valueWithArguments: {func2. x. y. h. 'y'}.
^ det2x2 valueWithArguments: {fx. dfy1. fy. dfy2}.
].
bN := [:func1 :func2 :x :y :h |
| dfx1 fx dfx2 fy |
dfx1 := partialDerivative valueWithArguments: {func1. x. y. h. 'x'}.
fx := func1 value: x value: y.
dfx2 := partialDerivative valueWithArguments: {func2. x. y. h. 'x'}.
fy := func2 value: x value: y.
^ det2x2 valueWithArguments: {dfx1. fx. dfx2. fy}.
].
" Функция для вычисления метода Ньютона "
newtonMethod := [:func1 :func2 :x0 :y0 :h :epsilon |
| x y continueIteration |
x := x0.
y := y0.
continueIteration := true.
[continueIteration] whileTrue: [
| jacobianValue aNValue bNValue newX newY |
jacobianValue := jacobian valueWithArguments: {func1. func2. x. y. h}.
aNValue := aN valueWithArguments: {func1. func2. x. y. h}.
bNValue := bN valueWithArguments: {func1. func2. x. y. h}.
newX := x - (aNValue / jacobianValue).
newY := y - (bNValue / jacobianValue).
((newX - x) abs <= epsilon and: [(newY - y) abs <= epsilon])
ifTrue: [continueIteration := false]
ifFalse: [x := newX. y := newY].
].
^ Array with: x with: y.
].
" Начальные значения и параметры "
x0 := 1.0.
y0 := 1.0.
h := 1e-6.
epsilon := 1e-10.
" Вычисление метода Ньютона "
| result |
result := newtonMethod valueWithArguments: {f. g. x0. y0. h. epsilon}.
result ifNotNil: [
Transcript show: 'x = ', result first printString, ', y = ', result second printString; cr.
] ifNil: [
Transcript show: 'Метод Ньютона не сошелся'; cr.
].