구 표면 위의 벡터 필드 그리기

먼저 1 C(쿨롱)과 -1 C의 전하가 각기 (0,0,-1)과 (0,0,1)의 위치에 있을 때, 이들이 생성하는 전기장을 구해보자. 이들 두 전하가 만드는 포텐셜을 구하고, 이 포텐셜의 음의 기울기를 구하면 전기장을 알 수 있다.

φ[x, y, z] == 1/(4 π Subscript[ε, 0])
(1/Norm[{x, y, z} - {0, 0, -1}] - 
  1/Norm[{x, y, z} - {0, 0, 1}])
% // Simplify[#, x_ \[Element] Reals] &
% // Map[-HoldForm@Grad[#, Cartesian[x, y, z]] &, #] & // 
  ReplaceAll[#, -HoldForm@
       Grad[φ[x, y, z], Cartesian[x, y, z]] -> 
     OverVector[Ε]] & // ReleaseHold
eq1 = %;

편의를 위해, 공통 계수를 없애면 아래와 같은 벡터 필드 방정식을 얻을 수 있다.

4 π Subscript[ε, 0] eq1[[2]]
eq2 = %;

이 수식을 VectorPlot3D 함수를 이용해서 필드를 그리면, 직육면체 영역 안에서 전기장을 나타내는 벡터들을 보여준다.

eq2
% // VectorPlot3D[#, {x, -1, 1}, {y, -1, 1}, {z, -1, 1},
   VectorScale -> {Automatic, Automatic, None},
   VectorColorFunction -> "Rainbow"] &

이제 원점에서 1의 거리에 있는 전기장만을 그리려면, 즉 반지름이 1인 원의 표면 위에 있는 벡터만을 그리려면 어떻게 해야 할까? ‘델타함수의 그래프 그리기’에서 한 것처럼 기본 Plot 계열의 함수가 처리해 주지 못할 때엔 ListPlot 계열의 함수를 이용하는 방법을 생각해본다.

Table[{r, θ, ф}, {r, 1, 1}, {θ, 0, π, π/
    10}, {ф, 0, 2 π, π/10}] // Flatten[#, 2] &;
% // Map[CoordinatesToCartesian[#, Spherical] &, #] &;
% // Map[{#, (4 π Subscript[ε, 0]
          eq1[[2]] /. {x -> #[[1]], y -> #[[2]],
         z -> #[[3]]})} &, #] &;
% // ListVectorFieldPlot3D[#, Axes -> True,
   AxesLabel -> {"x", "y", "z"}, VectorHeads -> True,
   ScaleFunction -> (0.2 &)] &

이 복잡한 명령어들이 하는 일은 다음과 같다. 우선, 첫 번째 줄에서 구 좌표계를 이용해서 구 표면 위의 일정한 좌표를 얻는다. 두 번째 줄에선 우리가 가진 식이 직교 좌표계의 식이므로 구 좌표계의 좌표를 직교 좌표계의 좌표로 바꾼다. 세 번째 줄에선 좌표를 식에 대입해서 그래프를 그릴 자료를 생성한다. 네 번째 줄에서는 ListVectorFieldPlot3D를 이용해서 전기장의 그래프를 그린다.

버전 7.0에선 ListVectorFieldPlot3D보다 ListVectorPlot3D를 사용하길 권하지만, ListVectorPlot3D에 아직 문제가 있는 듯이 보인다.

델타함수 그래프 그리기

코사인 함수의 푸리에 변환은 두 델타함수의 합으로 나타난다. 하지만, 그래프로 그리면 델타함수의 피크를 볼 수 없다.

FourierTransform[Cos[1.1 t], t, ω]
% // Plot[#, {ω, -3, 3}] &

코사인 함수의 푸리에 변환으로 얻은 함수에는 그래프로 그리기 어려운 두 가지 이유가 있다. 먼저, 델타함수는 인자를 0으로 만드는 지점에서 무한대로 발산하므로 그래프로 표현하기 어렵다. 관례적으로, 델타함수를 그래프로 그릴 땐 발산 점에서 델타함수의 값을 1로 대체해서 그래프를 그린다. 즉, DiracDelta를 DiscreteDelta로 치환하면 된다. 또한, 델타함수는 연속함수가 아니어서 Plot이 발산 점을 찾아내어 그려주지도 못하는데, 이 문제는 발산 점을 포함하는 목록을 만들어서 함수 자료에 발산 점이 꼭 들어가게 하여 해결한다.

FourierTransform[Cos[1.1 t], t, ω]
% /. DiracDelta[x_] :> DiscreteDelta[x]
% // Table[{ω, #}, {ω, -3, 3, 0.001}] & //
 ListPlot[#, Joined -> True] &

이번에도 그래프엔 아무것도 나타나지 않는다. 진법 변환 오차 때문에 발생하는 문제로, DiracDelta[1.1 - (1.099 + 0.011)]를 실행해보면 1이 아닌 0이 나오는 것을 알 수 있다. 컴퓨터는 1.100과 1.099+0.001을 다른 수로 본다. 십진수로 입력한 수는 일단 이진수로 바뀌는데, 십진수 실수 중에는 이진수로 변환했을 때 무한소수가 되는 수가 있다. 예를 들어, 1.100은 이진수로 표현하면 소수점 아래 두 번째부터 다섯 번째까지가 반복하는 무한소수가 된다.

BaseForm[1.100, 2]

컴퓨터는 저장공간의 제약으로 무한소수를 일정 자리에서 버림하여 저장하는데, 이때 입력받은 십진수와 저장된 이진수 사이에 차이가 발생하게 된다. Rationalize를 이용해서 모든 수를 정수로 표현하면 이러한 문제를 없앨 수 있다.

FourierTransform[Cos[1.1 t], t, ω]
% /. DiracDelta[x_] :> 
  Hold@DiscreteDelta@Rationalize@Unevaluated@x
% // Table[{x, # /. ω -> x}, {x, -3, 3, 0.001}] & //
  ReleaseHold // ListPlot[#, Joined -> True] &

그래프의 특정 부분을 강조하려면

그래프에서 특정 부분을 강조하려면 Prolog나 Epilog 옵션을 이용한다. 예를 들어 코사인 그래프에서 π부터 2π까지의 영역을 강조하려면 다음과 같이 Prolog 옵션을 준다.

Plot[Cos[x], {x, 0, 10}, 
 Prolog ->; {Green, Opacity[0.1], 
  Rectangle[{\[Pi], -10}, {2 \[Pi], 10}]}]

Epilog를 이용하면 그래프를 가리게 되는데, 이때엔 Opacity 명령어로 투명도를 조절하면 숨겨진 그래프를 볼 수 있다. 이 방법은 Graphics로 표현되는 모든 그래프에서 사용할 수 있다.

GIST의 위도와 경도

구글 지도가 제공하는 OpenAPI인 구글 Geocoding API를 이용하면, GIST의 위도와 경도를 찾을 수 있다. 웹브라우저의 주소 입력 칸에 아래 주소를 넣으면, XML 형태의 위치 정보를 받을 수 있다.

http://maps.google.com/maps/api/geocode
/xml?address=GIST,+Oryong-dong+1,+Buk-gu,+Gwangju&
sensor=false

아래 결과를 보면, geometry 태그 아래에 위도와 경도 정보가 나와 있는 것을 알 수 있다.

구글 Static Maps API를 이용해서 우리가 얻은 위도와 경도 정보가 정확한지 확인해 볼 수 있다.

http://maps.google.com/maps/api
/staticmap?center=35.2278792,126.8415491&zoom=15&
markers=35.2278792,126.8415491&size=400 x400&sensor=false

이상의 과정을 매스매티카에서도 쉽게 재현해볼 수 있다. 먼저, XML 패키지를 읽어들인다.

<< XML`

구글 지도에서 위치 정보를 XML로 받아온다. 매스매티카는 웹브라우저와 달리 언어 정보를 보내지 않으므로, 결과를 한글로 받아오려면 언어 정보를 명시해줘야 한다.

gistXML =
 XMLGet["http://maps.google.com/maps/api/geocode
/xml?address=GIST,+Oryong-dong+1,+Buk-gu,+Gwangju&
language=ko&sensor=false"]

받아온 자료에서 위도와 경도 정보를 빼낸다.

latlng = ToExpression@
  Flatten@Cases[gistXML,
    XMLElement[
      "location", _, {XMLElement["lat", _, {lat_}],
       XMLElement["lng", _, {lng_}]}] -> {lat, lng}, Infinity]

이렇게 얻은 위도와 경도 정보를 지도로 확인해 볼 수 있다.

Block[{url, loc},
 loc = ToString /@ latlng /. {x_, y_} :> x <> "," <> y;
 url = "http://maps.google.com/maps/api/staticmap?center=" <> loc <>
   "&zoom=15&markers=" <> loc <> "&size=400x400&sensor=false";
 Import[url]
 ]

참고자료

1. The Google Geocoding API
   
2. Static Maps API V2 Developer Guide

이맥스의 파이썬 모드에서 주석처리

이맥스의 파이썬 모드(M-x python-mode)를 사용하다 보면, 중요한 기능 하나가 메뉴에 없는 것을 알 수 있다. 바로 주석 처리 기능인데, 사실 메뉴에만 없지 이미 모든 프로그램 모드에서 제공하는 기능이다. 단축키는 M-; (‘coent-dwim’)이다.

일정 부분을 주석 처리하려면, 주석처리를 시작할 부분에서 C-SPC-SPC로 Transient-mark-mode를 켜고 끝 부분까지 커서를 이동한다. 다음으로, M-;을 누르면 해당 구역이 주석 처리된다. 이미 주석 처리된 부분에서 이렇게 하면 주석이 해제된다.

임의의 줄에서 M-;을 누르면 줄의 맨 끝에서 주석을 시작한다.

정적분 입력 서식의 단축키

매스매티카에서 정적분(definite integral)을 입력할 때, 단번에 입력 서식을 만들 수 있는 단축키를 제공한다. Esc+dintt+Esc를 입력하면, 다음과 같은 입력 서식이 생성된다.

파이썬 스크립트에 한글 넣기

파이썬 스크립트에 한글을 쓰려면 아래와 같이 인코딩 설정을 해 주어야 한다. 코드 부분뿐만 아니라 주석 부분의 한글 입력에도 인코딩 설정이 필요하다. 인코딩 설정은 인터프리터 지시어 바로 다음 줄에 와야 한다.

#!/usr/bin/env python
# -*- coding: utf-8 -*-
...

템플릿 프로그래밍에서 복소수를 실수 또는 복소수형 변수에 대입

템플릿으로 생성한 클래스에 적합한 변수형을 인자로 받아서 내부에 저장된 복소수 값을 이용하여 인자를 변경하는 상황을 가정해보자. 아래의 ClassA 클래스에 있는 update 메서드는 객체에 적합한 변수형을 인자로 받지만, 내부 계산은 복소수로 처리한다.

#include 

template 
class ClassA
{
public:
    ClassA(): cmplx(0,0) {}
    void update(T& out) { out = cmplx; }
private:
    std::complex cmplx;
};

문제는 복소수를 인자에 대입하는 부분에서 발생한다. out이 복소수형일 때엔 문제가 안 되지만, float이나 double과 같은 실수형이라면 형변환 에러가 발생한다. STL의 std::complex 타입은 실수형으로의 변환 연산을 제공하지 않기 때문이다. 즉, std::comlex<double> 형으로 생성한 ClassA는 별문제가 없지만,

ClassA<std::complex<double> > obj;
std::complex<double> var;

obj.update(var);

double 형으로 생성한 ClassA는

ClassA obj;
double var;

obj.update(var);

컴파일하면, 다음과 같은 에러가 발생한다.

error: cannot convert ‘std::complex<double>’ to ‘double’ in assignment

이 문제는 복소수형 변수를 실수형 변수로 변환해주는 함수를 만들어 해결한다. ClassA::cmplx의  실수부에 원하는 값이 저장되어 있다고 한다면, std::complex.real()의 반환값을 실수형 변수에 대입하면 된다.

template <typename S, typename T>
static inline T& assign(const std::complex<S>& in, T& out)
{
    return out = static_cast<T>(in.real());
}

하지만 이 함수를 사용하면, 오히려 복소수형 변수를 복소수형 변수에 대입하는 경우가 문제가 된다. 따라서 복소수형에 대해서는 다음과 같이 부분 템플릿 특수화(partial template specialization)를 이용하여 실수부와 허수부끼리 대입해주어야 한다.

template <typename S, typename T>
static inline std::complex<T>& assign(const 
std::complex<S>& in, std::complex<T>& out)
{
    out.real(static_cast<T>(in.real()));
    out.imag(static_cast<T>(in.imag()));
return out;
}

이제 update 함수는 아래와 같이 바꾸면 형변환 문제가 해결된다.

...
    void update(T& out) { assign(cmplx, out); }
...