Tìm bán kính hình cầu dựa trên phương pháp bình phương nhỏ nhất

[lehiep232]

Xin chào mọi người.
Em có 1 bài toán.
Sau khi thu thập số liệu từ máy toàn đạc, yêu cầu là tìm ra bán kính và toạ độ tâm hình cầu.
Cách làm như sau:
Giả sử có 4 điểm (thực tế khoảng 40 hoặc nhiều hơn) như trong hình
Ban đầu ta giả thuyết rằng tâm hình cầu có toạ độ là (a1,b1,c1) = (0,0,0) bán kính r1=0
B1: Ta sẽ tính được khoảng cách từ tâm đến các điểm theo công thức = sqrt{ (0-xi)^2 + (0-yi)^2 + (0-zi)^2 }
B2: Tính tổng các giá trị trên chia cho số điểm ta được bán kính mới là r2
B3: Tính các giá trị (xi-a)/ri, (yi-b)/ri, (zi-b)/ri
Tổng từng giá trị trên
B4: Tìm được tâm mới theo công thức a2 = sum(xi)-r2*sum((xi-a)/ri))/n, tương tự cho b2, c2

Sau khi có r2, và tâm mới (a2, b2. c2) lặp lại B1-B4 để có r3 và tâm mới (a3, b3, c3)
Lặp lại quá trình trên cho đến khi nào độ lệch của ri và ri+1 là 0.00001

Nhờ mọi người viết giúp em công thức vòng lặp để tính cho nhanh.
Em cảm ơn.

 

[HieuCD]

Xin chào mọi người.
Em có 1 bài toán.
Sau khi thu thập số liệu từ máy toàn đạc, yêu cầu là tìm ra bán kính và toạ độ tâm hình cầu.
Cách làm như sau:
Giả sử có 4 điểm (thực tế khoảng 40 hoặc nhiều hơn) như trong hình View attachment 276029
Ban đầu ta giả thuyết rằng tâm hình cầu có toạ độ là (a1,b1,c1) = (0,0,0) bán kính r1=0
B1: Ta sẽ tính được khoảng cách từ tâm đến các điểm theo công thức = sqrt{ (0-xi)^2 + (0-yi)^2 + (0-zi)^2 }
B2: Tính tổng các giá trị trên chia cho số điểm ta được bán kính mới là r2
B3: Tính các giá trị (xi-a)/ri, (yi-b)/ri, (zi-b)/ri
Tổng từng giá trị trên
B4: Tìm được tâm mới theo công thức a2 = sum(xi)-r2*sum((xi-a)/ri))/n, tương tự cho b2, c2

Sau khi có r2, và tâm mới (a2, b2. c2) lặp lại B1-B4 để có r3 và tâm mới (a3, b3, c3)
Lặp lại quá trình trên cho đến khi nào độ lệch của ri và ri+1 là 0.00001

Nhờ mọi người viết giúp em công thức vòng lặp để tính cho nhanh.
Em cảm ơn.

Dùng hàm tự tạo

Function HinhCau(ByVal Rng As Range) As Variant
  Dim arr(), aTotal#(1 To 3), t#(), e#, res#, tmp#
  Dim r#, a#, b#, c#, sR&, sC&, i&, j&

  e = 0.00001 'Sai so cho phep
  arr = Rng.Value
  sR = UBound(arr): sC = UBound(arr, 2)
  For j = 1 To sC
    For i = 1 To sR
      aTotal(j) = aTotal(j) + arr(i, j)
    Next i
  Next j
  Do
    tmp = res
    ReDim t(0 To sC)
    For i = 1 To sR
      r = Sqr((arr(i, 1) - a) ^ 2 + (arr(i, 2) - b) ^ 2 + (arr(i, 3) - c) ^ 2)
      t(0) = t(0) + r
      t(1) = t(1) + (arr(i, 1) - a) / r
      t(2) = t(2) + (arr(i, 2) - b) / r
      t(3) = t(3) + (arr(i, 3) - c) / r
    Next i
    res = t(0) / sR
    a = (aTotal(1) - res * t(1)) / sR
    b = (aTotal(2) - res * t(2)) / sR
    c = (aTotal(3) - res * t(3)) / sR
  Loop Until Abs(res - tmp) <= e
  HinhCau = Array(res, a, b, c)
End Function

Xem cách dùng hàm trong file

 

[lehiep232]

Dùng hàm tự tạo

Function HinhCau(ByVal Rng As Range) As Variant
  Dim arr(), aTotal#(1 To 3), t#(), e#, res#, tmp#
  Dim r#, a#, b#, c#, sR&, sC&, i&, j&

  e = 0.00001 'Sai so cho phep
  arr = Rng.Value
  sR = UBound(arr): sC = UBound(arr, 2)
  For j = 1 To sC
    For i = 1 To sR
      aTotal(j) = aTotal(j) + arr(i, j)
    Next i
  Next j
  Do
    tmp = res
    ReDim t(0 To sC)
    For i = 1 To sR
      r = Sqr((arr(i, 1) - a) ^ 2 + (arr(i, 2) - b) ^ 2 + (arr(i, 3) - c) ^ 2)
      t(0) = t(0) + r
      t(1) = t(1) + (arr(i, 1) - a) / r
      t(2) = t(2) + (arr(i, 2) - b) / r
      t(3) = t(3) + (arr(i, 3) - c) / r
    Next i
    res = t(0) / sR
    a = (aTotal(1) - res * t(1)) / sR
    b = (aTotal(2) - res * t(2)) / sR
    c = (aTotal(3) - res * t(3)) / sR
  Loop Until Abs(res - tmp) <= e
  HinhCau = Array(res, a, b, c)
End Function

Xem cách dùng hàm trong file

Anh xem giúp
Bị lỗi sao không chạy đúng nhỉ, kiểm tra kĩ rồi

 

[befaint]

Bị lỗi sao không chạy đúng nhỉ, kiểm tra kĩ rồi

Bạn cần bật cho phép macro chạy. Chưa biết bật thì tìm từ khóa: Bật macro Excel.

 

[lehiep232]

Mình đã test với 12 điểm, nhưng kết quả không đúng với kết quả tính thủ công.

 

[HieuCD]

Mình đã test với 12 điểm, nhưng kết quả không đúng với kết quả tính thủ công.

Kết quả của mình tính theo sai số cho phép là 0.00001 ứng với bước thứ 13 trong file
Nếu muốn chính xác hơn thì thay đổi thông số nầy trong code: e=0.00001

 

[lehiep232]

Nếu mình muốn:
+ Hoặc thêm ràng buộc là toạ độ tâm hình cầu (a, b, c) cũng có có sai số là 0.00001 thì code thay đổi sao nhỉ.
+ Hoặc bỏ ràng buộc là sai số 0.00001. Chỉ cần chạy 102 vòng cho ra kết quả
Mình cảm ơn.

 

[HieuCD]

Nếu mình muốn:
+ Hoặc thêm ràng buộc là toạ độ tâm hình cầu (a, b, c) cũng có có sai số là 0.00001 thì code thay đổi sao nhỉ.
+ Hoặc bỏ ràng buộc là sai số 0.00001. Chỉ cần chạy 102 vòng cho ra kết quả
Mình cảm ơn.

Số vòng chạy lệ thuộc vào đặc điểm của dữ liệu, file đầu chạy hơn 30.000 vòng mới ra kết quả nên con số cố định 102 vòng có rủi ro về sai số khá cao
Code xét sai số cho cả 4 tham số

Function HinhCau(ByVal Rng As Range) As Variant
  Dim a(), aTotal#(1 To 3), t#(), res#(0 To 3), b
  Dim r#, c#, sR&, sC&, i&, j&, e#

  e = 0.00001 'Sai so cho phep
  a = Rng.Value
  sR = UBound(a): sC = UBound(a, 2)
  For j = 1 To sC
    For i = 1 To sR
      aTotal(j) = aTotal(j) + a(i, j)
    Next i
  Next j
  Do
    b = res
    ReDim t(0 To sC)
    For i = 1 To sR
      r = Sqr((a(i, 1) - b(1)) ^ 2 + (a(i, 2) - b(2)) ^ 2 + (a(i, 3) - b(3)) ^ 2)
      t(0) = t(0) + r
      t(1) = t(1) + (a(i, 1) - b(1)) / r
      t(2) = t(2) + (a(i, 2) - b(2)) / r
      t(3) = t(3) + (a(i, 3) - b(3)) / r
    Next i
    res(0) = t(0) / sR
    res(1) = (aTotal(1) - res(0) * t(1)) / sR
    res(2) = (aTotal(2) - res(0) * t(2)) / sR
    res(3) = (aTotal(3) - res(0) * t(3)) / sR
  Loop Until Abs(res(0) - b(0)) <= e And Abs(res(1) - b(1)) <= e And Abs(res(2) - b(2)) <= e And Abs(res(3) - b(3)) <= e
  HinhCau = res
End Function

 

[VetMini]

Số vòng chạy lệ thuộc vào đặc điểm của dữ liệu, file đầu chạy hơn 30.000 vòng mới ra kết quả nên con số cố định 102 vòng có rủi ro về sai số khá cao
Code xét sai số cho cả 4 tham số
...

Nhiều tham số quá liệu có nhiều khả năng không hội tụ (convergence) ?

 

[HieuCD]

Nhiều tham số quá liệu có nhiều khả năng không hội tụ (convergence) ?

Đây là bài toán tính tổng bình phương khoảng cách các điểm đến 1 điểm cần tìm (tạm gọi là "tâm hình cầu") đạt giá trị nhỏ nhất, bài toán nầy dùng đạo hàm giải trực tiếp tìm hệ phương trình chuẩn tắc kết quả không có sai số, nhà toán học nào đó đã tìm cách giải trung gian bằng dãy số giá trị bán kính hội tụ, mình có cảm giác khi bán kính hội tụ và tiến dần đến "tâm hình cầu" thì các giá trị tọa độ tâm sẽ hội tụ

 

[PhanTuHuong]

Hóng kết quả bài toán thôi

 

[VetMini]

Đây là bài toán tính tổng bình phương khoảng cách các điểm đến 1 điểm cần tìm (tạm gọi là "tâm hình cầu") đạt giá trị nhỏ nhất, bài toán nầy dùng đạo hàm giải trực tiếp tìm hệ phương trình chuẩn tắc kết quả không có sai số, nhà toán học nào đó đã tìm cách giải trung gian bằng dãy số giá trị bán kính hội tụ, mình có cảm giác khi bán kính hội tụ và tiến dần đến "tâm hình cầu" thì các giá trị tọa độ tâm sẽ hội tụ

Tôi tin rằng nếu xét thêm tâm hình cầu sẽ làm tăng độ chính xác thì nhà toán học kia đã nói thẳng ra.
Việc tăng độ chính xác của bài toán hình cầu cũng có nhiều bài nghiên cứu trên các báo chí khoa học. Tôi tham khảo thử vài bài không thấy nói đến việc xét 4 con số như thớt đòi hỏi.

 

[HieuCD]

Tôi tin rằng nếu xét thêm tâm hình cầu sẽ làm tăng độ chính xác thì nhà toán học kia đã nói thẳng ra.
Việc tăng độ chính xác của bài toán hình cầu cũng có nhiều bài nghiên cứu trên các báo chí khoa học. Tôi tham khảo thử vài bài không thấy nói đến việc xét 4 con số như thớt đòi hỏi.

Bản thân thớt không rỏ về sai số nên đưa ra cách tính không chuẩn
"+ Hoặc thêm ràng buộc là toạ độ tâm hình cầu (a, b, c) cũng có có sai số là 0.00001 thì code thay đổi sao nhỉ.
+ Hoặc bỏ ràng buộc là sai số 0.00001. Chỉ cần chạy 102 vòng cho ra kết quả"
Cách tính thứ 2 càng tệ hơn