Giải hệ phương trình

[Ho Huu Tho]

Mình có một hệ hai phương trình như sau:
x = K(a-x)(b-x-y)
y =K'(b-y)(b-x-y)
Trong đó, x và y là ẩn, còn lại K, K', a, b là hằng số.
Bạn nào có kinh nghiệm có thể cho biết là về nguyên tắc thì có giải được hệ trên không?

 

[canh cut nd k16a]

Mình có một hệ hai phương trình như sau:
x = K(a-x)(b-x-y)
y =K'(b-y)(b-x-y)
Trong đó, x và y là ẩn, còn lại K, K', a, b là hằng số.
Bạn nào có kinh nghiệm có thể cho biết là về nguyên tắc thì có giải được hệ trên không?

Nếu

K, K', a, b là hằng số.

thì anh thay số cụ thể vào sẽ tính dễ hơn đấy ạ,
VD: K= 1, K' =2, a= 3, b=4
hệ: x= ( 3-x)(4-x-y)
y= 2(4-y)(4-x-y)
Giải hệ này quá đơn giản anh nhỉ

 

[lananh_vy_vp]

thế nếu ko thay số, canhcut có cách giải nào tổng quát ko?áp dụng đc cho tất cả các trường hợp ấy?

 

[pththao]

Mình có một hệ hai phương trình như sau:
x = K(a-x)(b-x-y)
y =K'(b-y)(b-x-y)
Trong đó, x và y là ẩn, còn lại K, K', a, b là hằng số.
Bạn nào có kinh nghiệm có thể cho biết là về nguyên tắc thì có giải được hệ trên không?

Về nguyên tắc phương trình trên tương đương với phương trình bậc ba nên giải được tổng quát theo công thức Cardano, hoặc dùng mathematica là nhanh nhất:
Đưa hệ về dạng:
t = b - x - y = x/(K(a-x)) = y/(K'(b-y)),
trong đó ta đặt biến mới t. Từ hai hệ thức cuối giải được x và y theo t,
x = aKt/(aKt+1), y = bK't/(bK't+1).
Thế x và y vào phương trình đầu,
t = b - aKt/(Kt+1) - bK't/(K't+1),
hay
a + t = b/(K't+1) + a/(Kt+1).
Phương trình trên sau khi qui đồng mẫu số sẽ là bậc ba. Cũng có thể có thủ thuật nào giải nhanh phương trình đó nhưng mình chưa nghĩ ra

 

[thaibeouu]

Mấy cái này cứ cố gắng nhờ Maple hay Mathematica nó giải cho các bác ạ

 

[Ho Huu Tho]

Mình chưa biết dùng mathematica, bạn nào chỉ giúp được không. Có phải vào trang web của họ để down phần mềm về và cài chương trình vào máy tính không?

 

[pththao]

Mathematica rất phổ biến chắc chắn tìm được trên torrent mình gửi hôm trước! Trong đó sẽ có cả keycrack cho student version. Bạn vào help của mathematica, xem Solve[], hoặc xem qua mạng:
http://reference.wolfram.com/mathematica/ref/Solve.html
Viết đúng theo ví dụ đó là ok thôi.

Phiên bản mã nguồn mở mình được giới thiệu là:
http://maxima.sourceforge.net/
http://www.sagemath.org/
http://www.gnu.org/software/octave/
nhưng chưa thử bao giờ.

PS: down ở Torrent tốt nhất là ta ra hàng để down cho chắc ăn

 

[thaibeouu]

Eo ơi em mê mẩn Mathematica :mrgreen
Download:
http://hotfile.com/dl/103855565/60ec025/Mathematjca.8.part1.rar.html
http://hotfile.com/dl/103856057/c6e8cb9/Mathematjca.8.part2.rar.html
http://hotfile.com/dl/103858190/e525997/Mathematjca.8.part3.rar.html
http://hotfile.com/dl/103858589/e7002c7/Mathematjca.8.part4.rar.html

hoặc:
http://www.fileserve.com/file/xvRfMTw
http://www.fileserve.com/file/Tp5N8pm
http://www.fileserve.com/file/zNc9pZV
http://www.fileserve.com/file/6Kbgaut

Yên tâm links em check cả r

 

[canh cut nd k16a]

Về nguyên tắc phương trình trên tương đương với phương trình bậc ba nên giải được tổng quát theo công thức Cardano, hoặc dùng mathematica là nhanh nhất:
Đưa hệ về dạng:
t = b - x - y = x/(K(a-x)) = y/(K'(b-y)),
trong đó ta đặt biến mới t. Từ hai hệ thức cuối giải được x và y theo t,
x = aKt/(aKt+1), y = bK't/(bK't+1).
Thế x và y vào phương trình đầu,
t = b - aKt/(Kt+1) - bK't/(K't+1),
hay
a + t = b/(K't+1) - a/(Kt+1).
Phương trình trên sau khi qui đồng mẫu số sẽ là bậc ba. Cũng có thể có thủ thuật nào giải nhanh phương trình đó nhưng mình chưa nghĩ ra

Em là hs PT nên em làm theo cánh thủ công thế này:
khỏi cần đặt t:
-Xét x=0 thì ko tìm được y => x khác 0
- Xét y=0 cũng chẳng tìm được x => y khác 0
=> (b-x-y) cũng khác 0 nốt
sau đó lấy tỉ lệ x/y =(K(a-x) ) / (K'(b-y)
tiếp đó là tìm ra mối qua hệ giữa x và y rồi thế vào một trong 2 phương trình ban đầu --> được pt bậc 3 --> giải
Bài này không nhất thiết phải đặt thêm ẩn làm gì cho rắc rối , anh (chị) pththao nhỉ

 

[pththao]

Em là hs PT nên em làm theo cánh thủ công thế này:
khỏi cần đặt t:
-Xét x=0 thì ko tìm được y => x khác 0
- Xét y=0 cũng chẳng tìm được x => y khác 0
=> (b-x-y) cũng khác 0 nốt
sau đó lấy tỉ lệ x/y =(K(a-x) ) / (K'(b-y)
tiếp đó là tìm ra mối qua hệ giữa x và y rồi thế vào một trong 2 phương trình ban đầu --> được pt bậc 3 --> giải
Bài này không nhất thiết phải đặt thêm ẩn làm gì cho rắc rối , anh (chị) pththao nhỉ

Hai cách đó là một mà, đặt ẩn để viết phương trình cho gọn và đẹp thôi, bản chất của t chỉ là kết nối x và y thôi ^^.

 

[Ho Huu Tho]

Sau khi thay số và nhập vào Mathematica hệ phương trình trên:
Solve[10^20*(10^-18 - x)*(10^-18 - x - y) == x &&
10^6 (10^-6 - y)*(10^-18 - x - y) == y && x > 0 && y > 0, {x, y}];
N[%]
Nghiệm của hệ mà chương trình đưa ra là:
{{x -> 1.*10^-18 + 7.42188*10^-34 I,
y -> 2.*10^-6 - 3.71094*10^-34 I}
{x -> 8.68286*10^-19 - 1.03762*10^-20 I,
y -> 6.58568*10^-20 + 5.18808*10^-21 I}}
Mình không hiểu chữ i trong thành phần của nghiệm là gì, bạn nào giải thích giúp được không?

 

[pththao]

Chữ I là số ảo. Thọ nên rescale biến trước khi giải số, chẳng hạn chia tất cả cho 10^-6 để tránh sai số nhỏ và lớn trong máy tính. Có thể đặt option cho biến và tham số là thực theo ví dụ.
Mathematica cũng như mọi algebra system không thể hoạt động một cách vạn năng. Nói chung nên tìm cách đơn giản hóa phương trình trước khi đưa vào. Cậu có thể thu được công thức tổng quát với các tham số a, b khi đơn giản hóa phương trình hợp lý. Chẳng hạn mình nghĩ sẽ có công thức tường minh cho t.

 

[pththao]

Code sau chạy trên R giải iteratively phương trình đó. Thọ có thể dùng để kiểm tra dịch chuyển cân bằng. Nếu cần điều chỉnh tham số hợp lý thì chỉ cẩn khéo léo đổi biến để không đưa vào các số lớn và quá bé thôi. Bạn save với đuôi mở rộng .R và chạy từ cmd của R bằng lệnh source('program.R') để chạy. Nếu xấp xỉ không tốt nó sẽ báo và tự thoát.

#HISTORY:
#Febuary 24, 2011: Created.
rm(list = ls())
#Function iteratively solve an eq.
#t0 = estimation of solution.
#a,K1,b,K2: parameters (rescaled).
#delta: error threshold.
itersolve = function(t0,a,b,K1,K2, delta = 0.0001)
{
N = 100000
t1 = t0
BAD = TRUE
for (k in 1:N)
{
t2 = a/(K1*t1+1) + b/(K2*t1+1) - a
if (abs(t2 - t1) > delta)
{
t1 = t2
}
else
{
BAD = FALSE
break
}
}
if (BAD)
{
stop('Solution was badly approximated, function exited')
}
return(t2)
}
#Parameters and first approximation:
t0 = 0
a = 2
b = 1
K1 = 0.1
K2 = 0.4
delta = 0.00001
#Solving:
t0 = itersolve(t0 ,a ,b , K1, K2, delta)
x = a*K1*t0/(K1*t0+1)
y = b*K2*t0/(K2*t0+1)

Đoạn code sau có thể khảo sát biến thiên của t0, bạn có thể dùng để kiểm tra điều kiện cân bằng.

allB = seq(1.2,3, len = 10)
allt0 = rep(0,10)
for (k in 1:length(allA))
{
b = allB[k]
allt0[k] = itersolve(t0 ,a ,b , K1, K2, delta)
}
dev.new()
plot(allB, allt0, col = 'red')
dev.copy2pdf(file = 'out.pdf') # Save file if necessary.

Với bộ số trên code cho nghiệm tốt (ổn định.) Nếu thay đổi các tham số quanh giá trị đó từng bước nhỏ mình nghĩ không vấn đề gì.

 

[Ho Huu Tho]

Sau khi save vào working directory thành file Thao_giai_he_phuong_trinh, mình ghi lệnh:
> source('Thao_giai_he_phuong_trinh.R')
... rồi enter, mà không thấy cái gì hiện ra cả. Hiện chưa biết sai ở khâu nào, nhờ Thảo chỉ giúp.

 

[pththao]

Ah, không hiện ra là ok mà, nếu hiện ra dòng 'Solution was badly approximated, function exited' thì mới là có vấn đề. Bạn gõ tiếp lệnh x (enter) để có giá trị của x; tương tự y (enter) để có giá trị của y. (Mình viết theo thuật toán nêu ở trên, thông qua giải phương trình cho t, phương pháp xấp xỉ liên tiếp, chưa kiểm tra lại, Thọ có thể thay lại vào phương trình xem có đúng không.)

 

[Ho Huu Tho]

Mình thấy phương trình chuẩn rồi, nhưng chưa biết cách để rescale. Tuy nhiên, mình cứ thay trực tiếp các số sau:
a=10^(-7)
b=10^(-6)
K1=10^20
K2=10^6
Thấy nó báo nghiệm là:
> x
[1] 5e-08
> y
[1] 9.999e-07
Không biết nghiệm này có tin cậy được không và giá trị của nghiệm của phải là: x = 5*e^(-8) không? Hay là x = 5*10^(-8)?

 

[pththao]

x = K(a-x)(b-x-y)
y =K'(b-y)(b-x-y)

Trước tiên 5e-8 là 5*10^(-8), qui ước viết dấu chấm động chứ không phải e^(-8). Khi nhập số a = 10^(-7) bạn cũng có thể nhập a = 1e-7, chương trình sẽ nhanh hơn chút, nhưng không đáng kể.

Cách rescale đơn giản nhất như sau: vì nói chung x và y mong đợi là có giá trị gần với a, b.
Bạn đặt x' = x/b, y' = y/b, (tương đương với đổi đơn vị cho x, và y) khi đó giá trị của x, y mong đợi là cỡ đơn vị. Khi đó bạn chuyển sang phương trình mới:

Thay vào phương trình trên được:
x' = Kb (a/b - x')(1 - x' -y')
y' = K'b(1 - x')(1 - x' - y')
Vậy có thể dùng biến mới b -> 1, a -> a/b, K -> Kb, K' -> K'b. Khi đó có thể thấy ngay các giá trị mang số hợp lý hơn (tuy nhiên nên kiểm tra đơn vị của K1 và K2 cẩn thận, mình thấy hắn lớn quá.)

Cuối cùng để xem nghiệm có hợp lý không đơn giản nhất là nhìn vào t0 (gọi bằng lệnh >t0), nếu bản thân giá trị của t0 không bé hơn delta thì bước đầu có thể tin được. Cẩn thận hơn ta thay ngược vào phương trình tính thử xem sao, có thể viết một đoạn code ở dưới để thử lại

 

[Ho Huu Tho]

Trước tiên 5e-8 là 5*10^(-8), qui ước viết dấu chấm động chứ không phải e^(-8). Khi nhập số a = 10^(-7) bạn cũng có thể nhập a = 1e-7, chương trình sẽ nhanh hơn chút, nhưng không đáng kể.

Cách rescale đơn giản nhất như sau: vì nói chung x và y mong đợi là có giá trị gần với a, b.
Bạn đặt x' = x/b, y' = y/b, (tương đương với đổi đơn vị cho x, và y) khi đó giá trị của x, y mong đợi là cỡ đơn vị. Khi đó bạn chuyển sang phương trình mới:

Thay vào phương trình trên được:
x' = Kb (a/b - x')(1 - x' -y')
y' = K'b(1 - x')(1 - x' - y')
Vậy có thể dùng biến mới b -> 1, a -> a/b, K -> Kb, K' -> K'b. Khi đó có thể thấy ngay các giá trị mang số hợp lý hơn (tuy nhiên nên kiểm tra đơn vị của K1 và K2 cẩn thận, mình thấy hắn lớn quá.)

Cuối cùng để xem nghiệm có hợp lý không đơn giản nhất là nhìn vào t0 (gọi bằng lệnh >t0), nếu bản thân giá trị của t0 không bé hơn delta thì bước đầu có thể tin được. Cẩn thận hơn ta thay ngược vào phương trình tính thử xem sao, có thể viết một đoạn code ở dưới để thử lại

Giá trị của K1 và K2 mình tính từ công thức:
DeltaG = RTln(K)
Trong đó:
DeltaG được dự đoán từ trình tự của ADN, có đơn vị là cal/mol
R là hằng số khí = 1.9858775 (cal K−1 mol−1)
T là nhiệt độ phản ứng (55+273=328), đơn vị là K
Còn (K) là hằng số cân bằng, mình nghĩ sẽ có đơn vị là mol−1
Chẳng hạn phức hợp PA' có deltaG=15 Kcal = 15.000 cal. Thay số sẽ được:
15.000=1.9858775*328*ln(K) -> ln(K) =23 -> K=1e10
Không biết tính toán như vậy thì có bị nhầm chỗ nào không nhỉ?

 

[pththao]

Đơn vị của hằng số cân bằng nói chung thay đổi trong các trường hợp khác nhau, trong trường hợp này có vẻ là M^-1. Cái này mình cũng hơi mơ màng khi học về cân bằng động học.
Mình xem qua wiki thì thấy có lẽ phải như sau: http://en.wikipedia.org/wiki/Equilibrium_constant
Nồng độ đo trong phản ứng phải đưa về đơn vị chuẩn của nồng độ khi đo G, ví dụ đo G ở 1 at thì áp suất phải đo bằng at; đo G ở 1 M thì nồng độ khi áp dụng phải đo ở M. Vì vậy \delta G mà Thọ nói sẽ phải đi kèm với điều kiện phản ứng xảy ra, ví dụ at, nồng độ chuẩn, ví dụ M, sau đó ta đổi tất cả về đơn vị chung đó.

 

[Ho Huu Tho]

Hệ gồm hai phương trình phản ứng:
A + A' = AA'; hệ số cân bằng là K1
p + A' = PA'; hệ số cân bằng là K2

Nồng độ ban đầu của A = A' là a; của P là b.
Lúc cân bằng có nồng độ của AA' là x; của PA' là y.

-> nồng độ lúc cân bằng của A là a-x, của A' là a-x-y, của P là p-y.

Viết lại hệ phương trình (hệ 1) một lần nữa cho hệ lúc cân bằng:
x = K1(a-x)(a-x-y)
y = K2(b-y)(a-x-y)

Rescale lại ẩn x và y như sau:
x' = x/a
y' = y/a

Thay x', y' vào hệ trên sẽ được hệ mới (hệ 2) là:
x' = aK1(1-x')(1-x'-y')
y' = aK2(b/a-y')(1-x'-y')

Dùng code của Thảo để giải hệ này với các tham số được thay đổi như sau:

K1 -> K1'= aK1
K2 -> K2'= aK2
a -> a' = 1
b -> b' = b/a

Thay số cụ thể:
(K1; K2; a; b) =
1.00E+06; 1.00E+12; 5.12E-14; 1.00E-06

(K1'; K2'; a'; b') =
5.12E-08; 5.12E-02; 1.00E+00; 1.95E+07

Code được thay số như sau:
#HISTORY:
#Febuary 24, 2011: Created.
rm(list = ls())
#Function iteratively solve an eq.
#t0 = estimation of solution.
#a,K1,b,K2: parameters (rescaled).
#delta: error threshold.
itersolve = function(t0,a,b,K1,K2, delta = 0.00001)
{
N = 100000
t1 = t0
BAD = TRUE
for (k in 1:N)
{
t2 = a/(K1*t1+1) + b/(K2*t1+1) - a
if (abs(t2 - t1) > delta)
{
t1 = t2
}
else
{
BAD = FALSE
break
}
}
if (BAD)
{
stop('Solution was badly approximated, function exited')
}
return(t2)
}
#Parameters and first approximation:
t0 = 0
a = 1
b = 1.95E+07
K1 = 5.12E-08
K2 = 5.12E-02
delta = 0.00001
#Solving:
t0 = itersolve(t0 ,a ,b , K1, K2, delta)
x = a*K1*t0/(K1*t0+1)
y = b*K2*t0/(K2*t0+1)

Code này cho nghiệm là:
> x
[1] 0.0009977034
> y
[1] 19480494
> t0
[1] 19505.86
Thử lại vào hệ 2 thấy không đúng. Nhìn qua nghiệm thì cũng thấy không phù hợp vì y' = y/a, nên nghiệm y' phải có giá trị nhỏ hơn 1 chứ nhỉ.
Nhờ Thảo xem lại giúp nhé.