2012-06-06 6 views
2

방정식의 풀이 방정식 시스템을 풀 필요가 있습니다. 여기에서 A는 n 행 x 행렬이고 B는 n 행 n 행렬입니다. 1 벡터입니다. 2D 포아송 문제에 대한 이산화 (discretization) 행렬이므로, 5 개의 대각선 만이 널이 될 수 없다는 것을 알고 있습니다. Lapack은이 특정 문제를 해결하는 함수를 제공하지 않지만 방정식의 밴드 매트릭스 시스템, 즉 DGBTRF (LU 분해) 및 DGBTRS를 해결하는 기능을 제공합니다. 이제 5 개의 대각선은 주 대각선, 첫 번째 대각선 위 및 아래의 주 대각선 및 두 개의 대각선 위와 아래 대각선 대각선 wrt 주 대각선. 밴드 스토리지에 관한 lapack 문서를 읽은 후에 A를 밴드 저장 포맷으로 저장하기 위해 (3 * m + 1) -by-n 매트릭스를 만들어야한다는 것을 배웠다.이 매트릭스 AB를 호출하자. 이제 질문 :방정식의 줄무늬 행렬 시스템을 해결하십시오.

1) dgbtrs와 dgbtrs_의 차이점은 무엇입니까? 인텔 MKL은 둘 다 제공하지만 이유를 이해할 수 없습니다.

2) dgbtrf는 밴드 스토리지 매트릭스가 어레이가되도록 요구합니다. 행이나 열로 AB를 선형화해야합니까?

3)이 두 함수를 호출하는 올바른 방법입니까?

int n, m; 
double *AB; 
/*... fill n, m, AB, with appropriate numbers */ 

int *pivots; 
int nrows = 3 * m + 1, info, rhs = 1; 
dgbtrf_(&n, &n, &m, &m, AB, &nrows, pivots, &info); 
char trans = 'N'; 
dgbtrs_(&trans, &n, &m, &m, &rhs, AB, &nrows, pivots, B, &n, &info); 

답변

0
  1. 는 또한 DGBTRSDGBTRS_ 제공합니다. 그것들은 당신이 신경 쓰지 말아야 할 포트란 행정부입니다. dgbtrs을 호출하십시오 (이유는 일부 아키텍처에서는 포트란 루틴 이름에 밑줄이 붙지 않고 다른 이름은 대문자 또는 소문자 일 수 있습니다). 인텔 MKL #define이 올바른 경우는 dgbtrs입니다.

  2. LAPACK 루틴은 열 주요 행렬 (즉, 포트란 스타일)을 필요로합니다. 다른 열 다음에 열을 저장하십시오. 너가 사용해야하는 끈으로 동거 된 저장은 열심히이지 않는다 : http://www.netlib.org/lapack/lug/node124.html.

  3. 나에게 좋을 것 같지만 작은 문제는 사전에 시도해보십시오 (항상 좋은 생각입니다). 또한 0이 아닌 info을 처리해야합니다 (오류보고 방법).

더 나은 스타일 대신 일반 intMKL_INT을 사용하는 것입니다,이 권리 유형에 대한 형식 정의입니다 (일부 아키텍처에 따라 달라질 수 있습니다).

dgbtrf을 호출하기 전에 pivots에 메모리를 할당해야합니다.

+0

약 # 2, LU 분해를 망치지 않습니까? 내 말은, AB가 행 주요임을 dgbtrf에게 알리는 옵션이 없다는 것입니다. 게다가 두 가지 경우의 주요 치수를 어떻게 설정해야합니까? – Patrik

+1

예, 당신은 A'x = B를 풀 수 있습니다. 이것에 대한 답을 정정하겠습니다. 당신이 펑키 한 일을하지 않는다면 (인자로 서브 매트릭스를 전달하는 것과 같은) 선행 차원의 경우, 보통 행 수와 같습니다. 다시 말하지만, LAPACK 문제를 해결하는 가장 좋은 방법은 손으로 해결할 수있는 작은 문제를 시도하는 것입니다. –

0

이것은 주제가 아닐 수도 있습니다. 그러나 푸 아송 방정식의 경우 FFT 기반 솔루션이 훨씬 빠릅니다. 잠재적 인 필드의 2D FFT를 - (k^2 + λ^2)로 나눈 다음 IFFT를 수행하십시오. 람다는 k = 0에 대한 발산을 피하기 위해 작은 수입니다. 5- 대각 방정식은 포아송 방정식의 대역 제한적 근사이며, 미분 연산자를 유한 차분으로 근사합니다.

http://en.wikipedia.org/wiki/Screened_Poisson_equation