e^t

我是风的歌

正在浏览 science 里的文章

工作是不是因为无聊?学习是不是因为无聊?介个不晓得里~~

话说卖萌可以买到妹纸吗?In urgent need ing~

有图有代码

从一个如下格式的文件中读取球的坐标和半径,然后画出来,使用了vtk库

n

x1 y1 z1 r1

x2 y2 z2 r2

xn yn zn rn

编译命令如下:

gcc example.cpp -I/usr/include/vtk-5.4 -L/usr/lib -lvtkCommon -lvtkViews -o example


#include <vtkSphereSource.h>
#include <vtkSmartPointer.h>
#include <vtkAppendPolyData.h>
#include <vtkCamera.h>
#include <vtkPolyDataMapper.h>
#include <vtkActor.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkProperty.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkTransform.h>
#include <vtkDepthSortPolyData.h>

vtkSmartPointer<vtkAppendPolyData> GenerateSpheres(int ntheta,int nphi)
{
    vtkSmartPointer<vtkAppendPolyData> appendData =
        vtkSmartPointer<vtkAppendPolyData>::New();

    ifstream sphereFile("spheres.3D");
    int nSphere;
    sphereFile >> nSphere;
    std::cout << "Number of spheres: " << nSphere << std::endl;
    double coordSphere[nSphere][4];
    for(int i = 0;i<nSphere;i++)
        for(int j = 0;j<4;j++)
            sphereFile >> coordSphere[i][j];

    for(int i = 0;i<nSphere;i++)
    {
        vtkSmartPointer<vtkSphereSource> sphereSource =
            vtkSmartPointer<vtkSphereSource>::New();
        sphereSource->SetThetaResolution(ntheta);
        sphereSource->SetPhiResolution(nphi);
        sphereSource->SetRadius(coordSphere[i][3]);

        sphereSource->SetCenter(coordSphere[i][0],
                coordSphere[i][1],
                coordSphere[i][2]);
        sphereSource->Update();
        appendData->AddInputConnection(sphereSource->GetOutputPort());
    }
    return appendData;
}

int main (int argc, char *argv[])
{
  int theta = 16;
  int phi = 16;

  // Generate a translucent sphere poly data set that partially overlaps:
  vtkSmartPointer<vtkAppendPolyData> translucentGeometry =
    GenerateSpheres(theta, phi);

  // generate a basic Mapper and Actor
  vtkSmartPointer<vtkPolyDataMapper> mapper =
    vtkSmartPointer<vtkPolyDataMapper>::New();
  mapper->SetInputConnection(translucentGeometry->GetOutputPort());

  vtkSmartPointer<vtkActor> actor =
    vtkSmartPointer<vtkActor>::New();
  actor->SetMapper(mapper);
  actor->GetProperty()->SetOpacity(1); //change to make translucent !!!
  actor->GetProperty()->SetColor(1, 0, 0);
  actor->RotateX(0); // put the objects in a position where it is easy to see
                       // different overlapping regions

  // Create the RenderWindow, Renderer and RenderWindowInteractor
  vtkSmartPointer<vtkRenderer> renderer =
    vtkSmartPointer<vtkRenderer>::New();
  vtkSmartPointer<vtkRenderWindow> renderWindow =
    vtkSmartPointer<vtkRenderWindow>::New();
  renderWindow->AddRenderer(renderer);
  vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor =
    vtkSmartPointer<vtkRenderWindowInteractor>::New();
  renderWindowInteractor->SetRenderWindow(renderWindow);

  // Add the actors to the renderer, set the background and size
  renderer->AddActor(actor);
  renderer->SetBackground(1, 1, 1);
  renderWindow->SetSize(600, 400);

  // Setup view geometry
  renderer->ResetCamera();
  renderer->GetActiveCamera()->Zoom(2.2); // so the object is larger

  int success = EXIT_SUCCESS;

  // Initialize interaction
  renderWindowInteractor->Initialize();

  // Start interaction
  renderWindowInteractor->Start();

  return success;
}

 

vtk的一些基础特性(http://apps.hi.baidu.com/share/detail/23071009):

VTk是在三维函数库OpenGL的基础上,采用面向对象的设计方法发展起来的。它有两种不同的方式:图像模型和可视化模型。

1、图像模型:支持3D几何数据绘制,3D体数据绘制,2D几何文字,图像绘制。

采用了3D图形系统简单易用的特点,同时也采用了图形用户接口的方法。整个图形模型表现3D图形系统的本质特点,主要有9类基本对象:渲染控制器,渲染窗口,渲染器,灯光,照相机,角色,特性,映射,变换。

1)渲染控制器定义与设备无关的坐标计算方法,并创建渲染窗口。

2)渲染窗口管理显示设备上的窗口,一个或多个绘制方法可在渲染窗口上创建一个场景,渲染窗口是用户图形界面,其中包括设置渲染窗口的大小,产生立体显示效果等方法。

3)渲染器是管理光源照相机和绘制对象等的位置、属性等,提供了世界坐标系,观察坐标系和现实坐标系之间的转换。

4)灯光可在场景中照亮绘制对象,可通过调用参数改变控制灯光的状态、照射角度,照射强度,颜色等,并支持点光源和平行光源。

5)照相机是定义观察者的位置,聚焦点和其他相关属性,参数可有调用者根据需要设置。

6)(相当于OpenGL的模型变换)角色代表渲染场景中的绘制对象实体,通过参数调节可以设置角色的位置,方向,渲染特性,引用,纹理映射等属性,并可对角色进行缩放。

7)属性是说明几何物体的一些特性,实现三维图形真实感。

8)映射指定了渲染数据和图形库中基本图元之间的关系,一个或多个角色可以使用相同的映射,有多个参数可对其进行控制。

9)变换是一个放置4×4变换矩阵的堆栈,可以进行各种操作。

2、可视化模型:利用数据流程模型,由过程对象和数据对象组成。过程对象包括可视化流程的模块及算法,数据对象包括数据处理及当数据在网络中流动时对数据进行的操作。过程对象又分为:源对象,过滤器,映射。源对象是可视化流程的起点,源对象包括从文件中读入及程序内部产生的数据,过滤器接收从源对象来的输入数据,处理数据及输出数据,映射指定了基本图元与数据之间的接口,接收过滤输出的数据,并把数据映射为基本图元。

 

References:

[1] http://www.vtk.org/Wiki/VTK/Examples/Cxx/Visualization/CorrectlyRenderingTranslucentGeometry

[2] The official vtk documentation

没事儿熟悉一下C语言,林火模型其实物理学家用来研究SOC(self-organized criticality),也就是自组织临界性的玩意儿。说到其使用价值,把一个森林简化到这种地步,值得怀疑啊。下面用了狗血的周期边界条件(敢问编过科学计算程序的人,谁没有啃过面包圈?!),初始条件是全部是森林(100*100的网格)。图一表示几个闪电打下了之后,森林鸭梨很大,纷纷dead;图二是图一中的强大火势快没了,因为森林已经七零八落了;图三说明森林恢复力惊人,继续和大火做起了永不停息的斗争……

下面是C语言代码,直接gcc就行了。上面的图我用Mathematica画的。差点忘了,秀一下森林人口变化(刚开始是一万):

森林人口

森林人口

顺便表示对wordpress排版无力吐槽……


#include <stdio.h>

#include <stdlib.h>

#define p_fire 0.001

#define p_grow 0.02

#define ISFIRE 10

#define ISNULL 0

#define ISTREE 1

double ran1()

{

    return 1.0*rand()/RAND_MAX;

}

int main(int argc, const char *argv[])

{

    const int NX = 100;

    const int NY = 100;

    int forest[NX][NY];

    int forest_n[NX][NY];

    int i,j;

    int ide,idw,jdn,jds,sum;

    int n;

    FILE *pFile;

    pFile = fopen("forestfire.dat","w");

    // initialization

    for(i=0;i

        for(j=0;j

            forest[i][j]=1;

            //if(ran1() < 0.818) forest[i][j] = 0;

            //else forest[i][j] = 1;

    // save the initial forest

    for(i=0;i

        for(j=0;j

            fprintf(pFile,"%2d",forest[i][j]);

        fprintf(pFile,"\n");

    }

    for(n=1;n

    {

        for(i=0;i

        {

            for(j=0;j

            {

                ide = i + 1 > NX - 1 ? 0 : i + 1;

                idw = i - 1 < 0 ? NX - 1 : i - 1;

                jdn = j + 1 > NY - 1 ? 0 : j + 1;

                jds = j - 1 < 0 ? NY - 1 : j - 1;

                sum = forest[i][jdn] + forest[i][jds] + forest[ide][j] + forest[idw][j] +

                    forest[ide][jdn] + forest[ide][jds] + forest[idw][jdn] + forest[idw][jds];

                if(forest[i][j] == ISTREE)

                    if(sum>10) // neighbour is on fire!

                        forest_n[i][j] = ISFIRE;

                    else // no neighbour is on fire, lightning sets in however

                        forest_n[i][j] = (ran1()                else if(forest[i][j] == ISNULL) // empty cell has a chance to grow new tree

                        forest_n[i][j] = (ran1()                else // with fire, the tree is dead

                        forest_n[i][j] = ISNULL;

            }

        }

        // update

        for(i=0;i

            for(j=0;j

                forest[i][j]=forest_n[i][j];

        // store data

        for(i=0;i

            for(j=0;j

                fprintf(pFile,"%2d ",forest[i][j]);

            fprintf(pFile,"\n");

        }

    }

    fclose(pFile);

    return 0;

}

孤立子

抢沙发

孤立子是数学和物理上最美丽的存在之一(待续)

Semi-cylinder Vorticity

抢沙发

A tutorial case of the Gerris CFD package.

A right heading flow with U=1m/s enters the domain from the left boundary. Upper and lower boundaries are frictionless solid walls. The right boundary is a simple outlet with zero pressure gradient.

Vorticity street is generated by the instability of the trailing vortices behind the half cylinder.

Princeton Ocean Model

抢沙发

Just one default case coming with the model. The domain is a portion of East U.S. Coast.

1. two-dimensional velocity field (not scaled by magnitude, just to show the directions)

2. three-dimensional velocity field (scaled by magnitude)

3. sigma coordinate

爱因斯坦的那道题~
以大神之名,恐怕大部分地球人都听说过~

1、英国人住红色房子
2、瑞典人养狗
3、丹麦人喝茶
4、绿色房子在白色房子左面
5、绿色房子主人喝咖啡
6、抽Pall Mall 香烟的人养鸟
7、黄色房子主人抽Dunhill 香烟
8、住在中间房子的人喝牛奶
9、 挪威人住第一间房
10、抽Blends香烟的人住在养猫的人隔壁
11、养马的人住抽Dunhill 香烟的人隔壁
12、抽Blue Master的人喝啤酒
13、德国人抽Prince香烟
14、挪威人住蓝色房子隔壁
15、抽Blends香烟的人有一个喝水的邻居
谁养鱼?

===============================================================

首先要做的是把信息重新排列下~
9、 挪威人住第一间房
14、挪威人住蓝色房子隔壁
4、绿色房子在白色房子左面

8、住在中间房子的人喝牛奶
5、绿色房子主人喝咖啡
1、英国人住红色房子
15、抽Blends香烟的人有一个喝水的邻居
6、抽Pall Mall 香烟的人养鸟
7、黄色房子主人抽Dunhill 香烟
10、抽Blends香烟的人住在养猫的人隔壁
11、养马的人住抽Dunhill 香烟的人隔壁
12、抽Blue Master的人喝啤酒
13、德国人抽Prince香烟

2、瑞典人养狗
3、丹麦人喝茶

===============================================================

第一突破口:房子位置(9+14+4)。推出排列顺序:
黄-蓝-绿-白-红
或者:
黄-蓝-红-绿-白

第一个推测不成立,因为(5)和(8)在这个排列下矛盾;所以房子排列顺序为第二个。插入可以插入的条件并从列表中删除它们(为方便香烟名称替换如下:
Dunhill-Du Blends-Bl; Pall Mall-PM; Blue Master-BM Prince-Pr)
这时我们的解答状态:

2、瑞典人养狗
3、丹麦人喝茶
6、抽PM 香烟的人养鸟
10、抽Bl香烟的人住在养猫的人隔壁
12、抽BM的人喝啤酒
13、德国人抽Pr香烟
15、抽Bl香烟的人有一个喝水的邻居

黄-蓝-红-绿-白
挪-    -英-          ——国
——-奶-咖      ——饮
Dh                     ——烟
—-马-               ——宠

===============================================================

第二突破口:信息(15)。Bl香烟不可能在白房子,因为邻居饮品是咖啡,所以只能在蓝房子,并且挪威人喝水。

2、瑞典人养狗
3、丹麦人喝茶
6、抽PM 香烟的人养鸟
10、抽Bl香烟的人住在养猫的人隔壁
12、抽BM的人喝啤酒
13、德国人抽Pr香烟

黄-蓝-红-绿-白
挪-    -英-         ——国
水-    -奶-咖    ——饮
Dh-Bl-             ——烟
—-马-             ——宠

===============================================================

这时饮+烟信息几乎满了,找限制最强的信息就是饮+烟的组合——信息12。然后是信息13、然后信息2。这道题至此已经迅速崩溃~

黄-蓝-红-绿-白
挪-丹-英-德-瑞——国
水-茶-奶-咖-啤——饮
Dh-Bl-PM-Pr-BM——烟
猫-马-鸟 狗——宠

题记:
所谓突破口自然是爱因斯坦设定的时候故意留下的,而“做出来就是2%聪明的人”也只是玩笑罢了。就算不是玩笑,2%这么大比例——解出来的成就感也是null哈~

OK… I confess I got “Einstein did it” from an unverified source (an old magazine). But this is an question I think even Einstein should be interested in [and solve it in a sec :) ], so…

Let me set up the scene:

We have a clock with a hour-hand and a minute-hand. And we know at certain moments we can exchange the position of these two hands and still get a valid time. The most obvious example is 12:00 o’clock when those two hands overlap. An opposite example is 6:00 o’clock — you won’t get 12:30 if you switch the hands, but get a configuration that a functional clock can never achieve!

Now the question is: during a day (0:00 am ~ 12:00 am), how many times we can make valid switches? DON’T WASTE A CHANCE FOR FUN! CLOSE THIS PAGE AND THINK ABOUT IT YOURSELF!

I once solved it several years ago using pure algebra (equations) at that time. Yesterday when I was about to sleep, this question came back to my mind and after a feeling of a flash in my mind, I got a new solution – a solution with the earliest, and most elegant branch of math – geometry.

Although I am not, but I can imagine a qualified mathematian do NOT dive into questions, but would stand aside for a while to analyze: such as making abstractions, before getting started.

===*===My Abstraction: this problem actually has nothing to do with time, because if the key is time, then this problem can be restated with a digital clock – but it can’t. So I would say this problem is a geometry problem (the angle is the key) with certain time-like constrains. Without constrains the configration space would be a Cartesian product of [0,360]*[0,360]. Any configration, making valid time or not, belongs to it. Configurations making a real time comprise a ‘subspace’. What about ‘switch’? Very simple, as commonly occurs, switch means “reflection”.===*===

For valid time: beginning from each o’clock, everytime the minute hand makes a full circle (360 degrees), the hour-hand moves 30 degree – resulting in 12 green lines.

Next reflect those lines by the red line (y=x), we get:

The solution is the number of crossings of green and blue lines: 144-1.

You should minus one, shouldn’t you?