Comsol with Matlab 应用实例(3)

mxio

一直有个想法，做一个随机圆柱体填充到立方体中，效果如图：

计算几何的算法会更好，无奈翻了下发现都是公式，不知道怎么去用。于是用了个偷懒的方法——试探。当然这种试探的方法，效率就很低了。
我的方法主要用到了：在底面找个随机位置，然后随机方向上放置圆柱体，然后判断是否和之前的圆柱体相交。
从而达到建立互不干涉的圆柱体。

1. function out = randcylinder
   
2. % This is a script for Comsol(4.x) with Matlab.
   
3. % It product random position cylinders inside a block.
   
4. %
   
5. % Copyright to mxio at Shanghai University.
   
6. %
   
7. %                                                        mxio
   
8. %                                                19,July,2013
   
9. import com.comsol.model.*
   
10. import com.comsol.model.util.*
    
11. model = ModelUtil.create('Model');
    
12. model.modelNode.create('mod1');
    
13. model.geom.create('geom1', 3);
    
14. 
    
15. blk_length = 1;
    
16. blk_size = [blk_length blk_length blk_length];
    
17. cyl_radius = 0.02 * blk_length;
    
18. cyl_length = blk_length + 4.5 * cyl_radius;       % about 3*sqrt(2)*cyl_radius
    
19. large_blk_position_z = -1.5 * cyl_radius;
    
20. blk_position = [0 0 0];
    
21. cy_num_total = 50;
    
22. cy_num = 2;
    
23. rng('shuffle'); % random number seed
    
24. 
    
25. model.geom('geom1').feature.create('cyl1', 'Cylinder');
    
26. model.geom('geom1').feature('cyl1').set('r', num2str(cyl_radius));
    
27. model.geom('geom1').feature('cyl1').set('h', num2str(cyl_length));
    
28. model.geom('geom1').feature('cyl1').set('pos', [rand rand large_blk_position_z]);
    
29. model.geom('geom1').feature('cyl1').set('ax2', [180*rand-90 360*rand]);
    
30. model.geom('geom1').run('cyl1');
    
31. 
    
32. 
    
33. cy_name_cell= {'cyl1'};
    
34. while cy_num<=cy_num_total
    
35.     model.geom('geom1').feature.create('uni1', 'Union');
    
36.     model.geom('geom1').feature('uni1').selection('input').set(cy_name_cell);
    
37.     model.geom('geom1').run('uni1');
    
38.     cy_name = ['cyl',num2str(cy_num)];
    
39.     model.geom('geom1').feature.create(cy_name, 'Cylinder');
    
40.     model.geom('geom1').feature(cy_name).set('r', num2str(cyl_radius));
    
41.     model.geom('geom1').feature(cy_name).set('h', num2str(cyl_length));
    
42.     model.geom('geom1').feature(cy_name).set('pos', [rand rand large_blk_position_z]);
    
43.     model.geom('geom1').feature(cy_name).set('ax2', [180*rand-90 360*rand]);
    
44.     model.geom('geom1').run(cy_name);
    
45.     model.geom('geom1').feature.create('int1', 'Intersection');
    
46.     tmp_name = [cy_name {'uni1'}];
    
47.     model.geom('geom1').feature('int1').selection('input').set(tmp_name);
    
48.     model.geom('geom1').run('int1');
    
49.     if isempty(mphshowerrors(model))
    
50.         model.geom('geom1').feature.remove(cy_name);
    
51.     else
    
52.         cy_name_cell{cy_num}= cy_name;
    
53.         cy_num = cy_num + 1;
    
54.     end
    
55.     model.geom('geom1').feature.remove('uni1');
    
56.     model.geom('geom1').feature.remove('int1');
    
57. end
    
58. 
    
59. model.geom('geom1').feature.create('blk1', 'Block');
    
60. model.geom('geom1').feature('blk1').set('size',blk_size);
    
61. model.geom('geom1').feature('blk1').set('pos', blk_position);
    
62. 
    
63. model.geom('geom1').feature.create('uni1', 'Union');
    
64. model.geom('geom1').feature('uni1').selection('input').set(cy_name_cell);
    
65. model.geom('geom1').run('uni1');
    
66. 
    
67. model.geom('geom1').feature.create('int1', 'Intersection');
    
68. tmp_name = {'uni1', 'blk1'};
    
69. model.geom('geom1').feature('int1').selection('input').set(tmp_name);
    
70. model.geom('geom1').run('int1');
    
71. 
    
72. model.geom('geom1').runAll;
    
73. mphgeom(model);
    
74. 
    
75. out = model;
    

复制代码

mxio

2013.7.19

rader001

NoEureka 发表于 2013-9-23 09:04
请问为什么我直接运行的时候，报错说Undefined function 'mphshowerrors'。我用的是comsol4.2 ...

4.2的不能用这个函数，需要修改。

showaround

楼主乃真汉子也！收下学习了。

atomical

楼主的方法简直了，另辟蹊径啊！

NoEureka

请问为什么我直接运行的时候，报错说Undefined function 'mphshowerrors'。我用的是comsol4.2

longwen36

楼主发的真是好贴啊，抛砖引玉

SMT_Captain

学习一下！

不会生气的蜗牛

楼主真棒，真是长见识

plwhu

谢谢了。。。。。

mr_xumao

Thanks for your share！！！

古襄月

不错。。

Joice晗

感觉真的高大上·~~

SComsol

试一下:)

cuiziyang

优秀的人一直都在，像大佬学习！

cuiziyang

优秀，像大佬学习66666

pcstudying

优秀的人，不仅自己优秀，还会分享，让我们学习