| 377 | } |
| 378 | |
| 379 | double meanOfDiamond( |
| 380 | const uint x, |
| 381 | const uint y, |
| 382 | const uint halfStepLength, |
| 383 | const std::vector<std::vector<double> > &R) { |
| 384 | |
| 385 | double average; |
| 386 | |
| 387 | if (x == 0) { // At top edge of system |
| 388 | if (PBC) { |
| 389 | average = 0.25*( |
| 390 | R[x][y + halfStepLength] + |
| 391 | R[x][y - halfStepLength] + |
| 392 | R[x + halfStepLength][y] + |
| 393 | R[R.size() - 1 - halfStepLength][y]); |
| 394 | } |
| 395 | else { |
| 396 | average = (1.0 / 3.0)*( |
| 397 | R[x][y + halfStepLength] + |
| 398 | R[x][y - halfStepLength] + |
| 399 | R[x + halfStepLength][y]); |
| 400 | } |
| 401 | } |
| 402 | else if (y == 0) { // At left edge of system |
| 403 | if (PBC) { |
| 404 | average = 0.25*( |
| 405 | R[x][y + halfStepLength] + |
| 406 | R[x][R[0].size() - 1 - halfStepLength] + |
| 407 | R[x + halfStepLength][y] + |
| 408 | R[x - halfStepLength][y]); |
| 409 | } |
| 410 | else { |
| 411 | average = (1.0 / 3.0)*( |
| 412 | R[x][y + halfStepLength] + |
| 413 | R[x + halfStepLength][y] + |
| 414 | R[x - halfStepLength][y]); |
| 415 | } |
| 416 | } |
| 417 | else { |
| 418 | average = 0.25*( |
| 419 | R[x][y + halfStepLength] + |
| 420 | R[x][y - halfStepLength] + |
| 421 | R[x + halfStepLength][y] + |
| 422 | R[x - halfStepLength][y]); |
| 423 | } |
| 424 | |
| 425 | return average; |
| 426 | } |
| 427 | |
| 428 | double nonPBCbottomEdgeDiamonds(const uint x, const uint y, const uint halfStepLength, std::vector<std::vector<double> >& R) { |
| 429 | |